AD- A 148  5 


CD 

00 


FINAL  TECHNICAL  REPORT 
1  February  1981  -  31  May  1984 
DETAILED  OCEANIC  CRUSTAL  MODELING 
Contract  N  • .  N00014-76-C-1070 


CALIFORNIA  INSTITUTE  OF  TECHNOLOGY 


>,v 


FINAL  TECHNICAL  REPORT 
1  February  1981  -  31  May  1984 


Name  of  Contractor 
Contract  No: 

Effective  Date  of  Contract: 
Contract  Expiration  Date: 
Total  Amount  of  Contract: 
Title  of  Work: 

Program  Director: 


Principal  Investigator: 


California  Institute  of  Technology 

N00014-76-C-1070 

February  1,  1981 

May  31,  1984 

$240,000 

Detailed  Oceanic  Crustal  Modeling 
J.  M.  McKisic 

Program  Director  of  Ocean  Acoustics  Program 
Scientific  Office,  Code  4250  A 
Office  of  Naval  Research 
Arlington,  VA  22217 

(202)  696-4204 

Donald  V,  Helmberger 
(818)  356-6998 


This  docusMat  Bos  baa 

public  Nkai  omd  sol*;  its 
distibutioa  k  mUmited. 


DT»^ 

EU 

OEC  1  4  1984 


Seismological  Laboratory 
Division  of  Geological  and  Planetary  Sciences 
California  Institute  of  Technology 
Pasadena,  CA  91125 


TABLE  OF  CONTENTS 


Summary . . . ^ 


Long  Period  Wave  Propagation  in  Laterally 
Varying  Structure  , 


Upper  Mantle  Shear  Structure  Beneath  the 
Northwest  Atlantic  Ocean 


Applications  of  the  Transmitted  Kirchhoff- 
Helmholtz  Method  to  Transmitted  Bodv  Waves 
and  Possible  Structural  Effects  at  NTS 


Accession  For 

NT IS  GRAAI 
DTIC  TAB 

Unannounced 
l^tLf  lcation. 


i'L 


B'  _ _ _ 

\  Distrlbatlan/ 

J  Availability  Codas 
jAvall  aad/or 
Dlst  Speelal 


Vo-  -V-  «•  .-.v.v.v.v'  -  w  V  W'j'V’.-V .M' Aiv.v.vvivivy 


w!.*' 


1 


I.  Summary 

The  research  performed  under  this  contract  can  be  divided  into  3  main 
topics:  changes  in  existing  methods,  Cagniard  de-Hoop  and  WKBJ,  which  enable 
construction  of  synthetics  for  mixed  path  situations:  use  of  Long  period  SH  waves 
with  source  in  the  Northwest  Atlantic  and  receivers  on  the  northeast  coast  of 
North  America  to  derive  an  oceanic  upper  mantle  shear  velocity  model;  and  a 
technique  based  on  evaluating  the  Kirchoff-Helmholtz  integral  for  predicting  the 
effect  of  near  source  or  near  receiver  structure  complexity  on  far  field  p  waves. 

In  Section  II  we  assess  the  fact  that  recent  models  of  upper  mantle  struc¬ 
ture  based  on  long  period  body  waves  (WWSSN)  suggest  large  horizontal  gra¬ 
dients,  especially  in  shear  velocities.  Some  changes  in  existing  methods  are 
required  to  construct  synthetics  for  mixed  path  situations.  This  is  accomplished 
by  allowing  locally  dipping  structure  and  making  some  modifications  to  general¬ 
ized  ray  theory.  Local  ray  parameters  are  expressed  in  terms  of  a  global  refer¬ 
ence  which  allows  a  de-Hoop  contour  to  be  constructed  for  each  generalized  ray 
with  the  usual  application  of  the  Cagniard  de-Hoop  technique.  Several  useful 
approximations  of  ray  expansions  and  WKBJ  theory  are  presented.  Comparisons 
of  the  synthetics  produced  by  these  two  basic  techniques  with  known  solutions 
demonstrates  their  reliability  and  limitations. 

In  Section  III,  we  have  modeled  the  SK  motion  from  earthquakes  in  the 
northwest  Atlantic  ocean  to  derive  an  oceanic  upper  mantle  shear  velocity 
model.  The  signals  were  recorded  on  long-period  WWSSN  and  Canadian  network 
stations  on  the  east  coast  of  North  America.  This  data  indicates  a  fast  (4.75 
km/sec)  lid  of  about  100  km  thickness  in  the  older  western  Atlantic.  Given  the 
lid  structure,  the  waveforms  and  traveltimes  from  the  more  distant  data  put 
tight  constraints  on  the  shear  velocities  at  greater  depths.  The  velocity  below 
200  km  was  found  to  be  indistinguishable  from  a  model  of  the  East  Pacific  Rise 
(Grand  and  Kelmberger,  1983)  found  using  the  same  technique.  We  find  the 
Canadian  shield  to  be  faster  than  both  the  old  northwestern  Atlantic  and  the 
young  East  Pacific  Rise  to  about  400  km  depth.  No  variations  below  400  km  are 
necessary  to  explain  the  data. 

In  Section  IV,  we  extend  the  Kirchhoff-Helmholtz  integral  method  to  calcu¬ 
late  acoustic  potentials  which  transmit  through  three  dimensional  warped  boun¬ 
daries.  We  specify  the  potentials  on  an  arbitrary  surface  with  Snell's  law  and 
plane-wave  transmission  coefficients  and  numerically  integrate  their  contribu¬ 
tions  at  a  receiver  via  the  scalar  integral  representation  theorem.  The  method 
is  appropriate  for  modeling  precritical  transmitted  potentials.  Results  from  test 
models  compare  well  with  optical  solutions  for  transmissions  through  a  flat 
interface. 


ABSTRACT 


Recent  models  of  upper  mantle  structure  based  on  long 
period  body  waves  (WWSSN)  suggest  large  horizontal  gradients, 
especially  in  shear  velocities.  Some  changes  in  existing  methods 
are  required  to  construct  synthetics  for  mixed  path  situations. 

This  is  accomplished  by  allowing  locally  dipping  structure  and 
making  some  modifications  to  generalized  ray  theory.  Local  ray 
parameters  are  expressed  in  terms  of  a  global  reference  which 
allows  a  deHoop  contour  to  be  constructed  for  each  generalized 
ray  with  the  usual  application  of  the  Cagniard-deHoop  technique. 
Several  useful  approximations  of  ray  expansions  and  WKBJ  theory 
are  presented.  Comparisons  of  the  synthetics  produced  by  these 
two  basic  techniques  with  known  solutions  demonstrates  their  reli- 

a 

ability  and  limitations. 

1.  Introduction 

Considerable  progress  has  been  made  recently  in  speeding-up  the  syn¬ 
thesizing  of  seismograms  with  the  introduction  of  WKBJ  and  Gaussian  beam 
methods,  see  Chapman  (1976)  and  Cerveny  et  al.  (1983).  These  methods  have 
proven  highly  useful  in  generalizations  to  laterally  varying  structure,  especially 
at  high  frequency,  see  for  example  Frazer  and  Phinney  (i960).  However,  in  the 
construction  of  longer  periods  (long  period  WWSSN  seismograms)  we  many  times 
are  interested  in  more  complete  solutions,  since  the  beginning  portion  of  sur¬ 
face  waves  become  important,  see  Grand  and  Helmberger  (1984b).  A  complete 
set  of  ray  parameter  contributions  is  required  to  construct  seismograms  in  this 
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situation.  In  particular,  one  needs  to  consider  ray  paths  leaving  the  source  hor¬ 
izontally.  a  case  where  the  WKBJ  method  breaks  down.  We  can  avoid  this  prob¬ 
lem  by  applying  a  mixture  of  generalized  ray  theory,  GRT,  and  WKBJ  or  Disk  rays 
as  defined  by  Wiggins  (1976). 

A  simple  example  of  this  procedure  is  given  in  Fig.  1  where  we  show 
schematically  how  to  construct  the  step  response  for  a  smooth  velocity  model 
approximated  by  a  stack  of  homogeneous  layers.  We  suppose  that  a  velocity 
model  can  be  chosen  such  that  the  step  response  remains  a  step  at  all  receiver 
positions.  The  simulation  of  this  step  can  be  achieved  by  summing  the  response 
from  three  energy  paths;  namely,  the  direct,  the  reflected  from  just  below  the 
source  or  reference  plane,  and  the  diving  WKBJ  contribution.  All  three  paths 
contain  a  product  of  the  transmission  coefficients  above  the  source.  The  WKBJ 
path  includes  the  transmission  coefficients  across  the  reference  plane,  taken  as 
the  interface  below  the  source.  We  have  included  a  diagram  of  the  0(f)  vs.  p 
curve  in  Fig.  1  for  reference  as  it  clearly  shows  that  the  diving  path  contributes 
little  except  at  the  larger  distances.  At  the  nearest  distance,  position  1,  the 
direct  ray  dominates.  The  reflected  path  contributes  some  as  critical  angle  is 
approached.  At  still  larger  distances,  position  3,  a  head  wave  along  the  bottom 
of  the  reference  interface  develops  followed  by  the  critically  reflected  pulse. 
Note  that  the  WKBJ  path  turns  on  at  the  start  of  the  headwave.  At  large  ranges, 
the  WKBJ  path  becomes  increasingly  dominant.  Note  the  interesting  interaction 
of  the  reflected  plus  the  direct  rays  where  the  high  frequency  energy  turns 
around  first.  This  frequency  dependent  behavior  forces  the  use  of  Airy  functions 
or  some  higher  order  expansion  near  the  turning  point.  In  practice,  this  step 
simulated  by  summing  these  three  energy  paths  does  not  quite  generate  a  step 
as  we  discuss  next. 


In  testing  the  accuracy  of  the  above  procedure  it  is  quite  useful  to  generate 
the  step  response  for  models  for  which  the  answer  is  known.  Thus,  we  begin  with 
a  homogeneous  fluid  whole  space  with  a  point  source  excitation  yielding  a  step 
response  at  all  positions  with  l/(distance)  decay.  We  next  impose  a  spherical 
coordinate  system  with  many  thin  shells  of  constant  velocity.  Applying  the  clas¬ 
sical  earth-flattening  approximation  we  obtain  a  model  with  a  smooth  velocity 
increase  in  depth,  see  Helmberger  (1973).  The  synthetics  generated  in  Fig.  2 
are  from  such  a  model  with  the  exact  step  responses  indicated  by  the  dotted 
lines  in  the  bottom  panel.  This  panel  also  displays  the  response  after  summing 
the  complete  set  of  generalized  rays;  direct  plus  rays  reflected  upward  from  all 
the  interfaces  below  the  source.  Chapman  (1976)  showed  that  the  sum  of  these 
rays  should  approach  a  step  for  the  diving  energy  portion  of  the  response  as  the 
layer  thicknesses  approach  zero.  However,  the  amplitude  will  be  (n/  3)  larger 
than  the  exact  answer  because  of  the  neglect  of  multiple  ray  interactions  with 
the  discrete  layering.  The  CRT  response  at  the  largest  distance  shows  the  most 
roughness  for  times  near  the  direct  arrival  when  the  interaction  with  the 
reflection  from  just  below  the  source  is  the  most  severe.  Similar  complexity 
occurs  with  the  hybrid  method  except  that  the  diving  energy  is  smoother  with 
WKBJ.  Short  period  synthetics  generated  from  these  step  responses  become 
quite  dirty  and  simple  geometric  ray  theory  yields  cleaner  results.  However,  for 
long  period  studies  the  advantage  of  being  able  to  include  the  radiation  pattern 
appropriate  for  earthquake  sources,  or  shear  dislocations,  far  outweighs  the 
disadvantage  of  the  noise  generated  by  the  hybrid  method.  For  example,  con¬ 
sider  the  SH-radiation  from  a  dip-slip  event  where  the  up-going  radiation  has 
opposite  polarity  from  the  down-going  energy,  see  Helmberger  (1973).  In  short, 
the  sum  as  displayed  in  Fig.  2  becomes  more  interesting  when  the  direct  ray 
trace  has  opposite  sign  from  the  other  two. 
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We  could  probably  improve  the  response  at  the  time  the  three  energy  paths 
interface  the  most  vigorously  by  including  a  few  more  GR’s  and/or  by  lowering 
the  reference  boundary  for  the  WKBJ  contribution.  However,  we  are  particularly 
interested  in  more  realistic  earth  models  with  a  crust-over-mantle  structure 
where  the  moho  or  low  velocity  zone  serve  as  natural  reference  boundaries. 
Thus,  we  propose  using  GR’s  to  compute  the  start  of  the  Love  waves,  and  WKBJ  to 
generate  the  responses  returning  from  deeper  structure.  This  approach  proved 
effective  in  studying  the  structure  and  evolution  of  the  lithosphere  for  an  old 
oceanic  plate,  Grand  and  Helmberger  (1984b).  It  would  be  advantageous  to 
treat  the  obvious  lateral  variation  encountered  in  such  studies.  Although  the 
real  world  is  truly  three  dimensional  in  nature,  some  useful  progress  can  be 
made  by  examining  profiles  of  data  along  paths  of  symmetry  where  two- 
dimensional  idealizations  are  appropriate.  We  will  address  such  models  in  this 
paper. 

Our  strategy  is  similar  to  Wiggins  (1976)  in  that  we  will  use  a  combination  of 
GRT  and  DRT  to  generate  synthetics  and  justify  the  latter  by  demonstrated  accu¬ 
racy. 

Review  of  Ray  interactions  with  nonplanar  structure 

Boundary  value  problems  involving  complicated  geometry  have  a  long 
rather  unrewarding  history;  thus,  we  will  jump  directly  to  approximate  solutions 
and  test  their  validity  against  finite-difference  calculations  and  other  more  well 
known  results.  Before  addressing  the  dipping  layer  problem  it  is  instructive  to 
examine  the  flat-layered  case  and  emphasize  the  geometric  interpretation  of 
generalized  ray  theory.  This  proves  particularly  useful  for  constructing  general¬ 
izations  to  more  complicated  situations  since  the  most  progress  in  understand¬ 
ing  these  problems  is  at  high  frequency.  Both  line  and  point  sources  will  be 


The  solution  of  the  scalar  wave  equation  assuming  line  source  excitation  for 
generalized  rays  as  given  by  Gilbert  and  Knopoff  (1961)  is 

=  H{t-t0)/  {tz-t?)wz  (1) 

where  t0  =  R/  a.R 2  =  za  +  ra,  and  a  =  velocity. 

$1  is  defined  as  the  displacement  potential  with  the  index  L  used  to  remind 
ourselves  of  the  line  excitation.  A  high  frequency  approximation  of  (l)  is 

aH(t-tg)/(t-t0)I/1!VW  (2) 

and  the  motion  decays  with  distance  as  the  v77.  The  solution  to  the  interface 
problem  setup  displayed  in  Fig.  3a  is 

$l(r,z,t)  =  lm(^—&-r(p))  (3) 

where 

f  =p(d,  +  dg)  +  At77,  +  hzrjs  (4) 

Vi  M-V-p*)1'* 

T  (p)  =  Transmission  coefficient 

The  symbol  (Im)  indicates  the  imaginary  part  of  the  complex  product  of  the 
functions  of  ray  parameter,  see  Helmberger  (1983)  for  example.  The  ray  param¬ 
eter  appropriate  for  the  direct  arrival  path,  p„,  can  be  obtained  by 
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we  see  that  the  ray  goes  from  the  source  to  the  receiver.  And 


sin0,  cos0j 
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t0  =  R\/  (Xj  +  R 2/  a2 


For  times  greater  than  f0,  we  must  solve  t  for  complex  p  such  that  the  ima¬ 
ginary  parts  of  pdx  and  77^  cancel,  etc. 

The  behavior  nearp0  can  be  approximated  by  noting  that 


t  =  +  -P.)  +  -p»)2/ 2 


and  solving  for 


(p  -pjz=  2(f  -f0)/(ff) 


Thus, 


Note  that  from  (4) 
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It  is  convenient  to  condense  the  various  factors  containing  p0  into 


Si(p.)» 


which  we  call  the  spreading  factor.  Thus, 


„  _  ^1  1/2 

i>L  =  - + - 5 - 

1hai  VzO-z[niaz\ 


¥e  note  that  by  differentiating  Snell's  law  we  obtain 


cos0,  COS0g 

- Lrf©  -  - 1^0 

ax  a2 


Substituting  this  expression  into  5^  we  obtain 


y/d&i 

R  xd0  i  +  R2d  0g 


cos0j 


If  Rz  -  0,  we  obtain  the  whole  space  spreading  again  where  (/Jjd©!)  is  just  the 
width  of  the  ray  tube  described  in  Fig.  3a.  A  correction  for  the  change  in  direc¬ 
tion  is  required  as  the  tube  crosses  the  interface,  namely 


(cOS0j/  COS0g)  . 


Thus,  the  denominator  of  (12)  is  again  the  width  of  the  ray  tube  at  the  receiver, 
L$  in  Fig-  3a.  Substituting  into  (3)  we  obtain 


*L  «  SLH(t  -t0)Re  (T(p0))/  (t  -t0)1/z 


where  Rt  indicates  the  real  part  operator. 


The  point  source  solution  for  the  same  pro  c'm  setup,  Fig.  3b,  is 


*P(r.2.t)  =  Vr£ 
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Vp  dp 
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and  applying  the  same  first  motion  approximation  we  obtain  a  slightly 
complicated  spreading  factor  namely. 


^sV  fv 

v  r  77, 


1  ,  daf 


1/2 


/F^in©,  +  A’asinBg 


sin©i 


-1/2 


Si 


and  note  that  by  letting  R2  =  0  we  obtain 


SP  =  1/  /? 


In  terms  of  area,  we  note  that 


SP  = 


sin©id  ©id$ 

(^isinQid^  +  tfasinQgdS’X/?,:*©,  +  /?ad©a 


cos©, 

COS0Z 


1/2 


(14) 

more 

(15) 

(16) 

(17) 

(IB) 


which  can  be  interpreted  as  the  incremental  element  of  area  at  the  source 
divided  by  the  projected  area  at  the  receiver  or  simply, 


Thus,  the  first  motion  behavior  becomes 


*p(rltlt)=  SpRe(T(p0))  (20) 

=  SpH(t-t0)Re(T(p0)) 

More  complicated  solutions  to  multi-layered  models  in  terms  of  ray  summations 
will  be  discussed  later. 


Although  CRT  for  parallel  interfaces  has  been  well  developed  the 
modifications  for  nonplanar  structure  or  smoothly  varying  interfaces  has  not. 
Some  of  the  difficulties  encountered  for  the  simple  wedge  problem  have  been 
discussed  by  Hudson  (1963).  Hong  and  Helmberger  (1977)  constructed  a  solu¬ 
tion  in  terms  of  generalized  rays  for  this  problem  and  defined  a  method  of  ray 
path  construction  compatible  with  the  usual  Cagniard-de  Hoop  formalism.  We 
will  consider  the  direct  arrival  interacting  with  two  dipping  interfaces  as  an 
example  application.  The  problem  setup  is  displayed  in  Figure  4a  with  the 
response  given  by 


*L  =Im  7’Ig(pl)TM(pz)^-  ^r 

r)i  at 


where  p  j  and  pj  are  defined  by  the  local  ray  parameter,  namely 


sin0. 

Pi  =  —  ,  p2  = 

f*i 


sin02 

<*2 


and  are  no  longer  equal.  However, 


- 11  - 


where 

02  =  02  +  ©J 

with  0,  defining  the  change  of  the  slope  of  interface  (l)  relative  to  the  previous 
reference  at  0j.  Performing  the  derivatives  discussed  in  the  previous  section  we 
obtain 


SL(Po)  = 


_ _ 

cos0,  cosQ2  cos©, 

R1dBl  +  R2dQ2 - —  +  ^3^03 - - - -7- 

COS02  COS02  COS0g 


which  is  similar  to  (12)  and  has  the  same  interpretation.  The  travel  time  is 
defined  by 


t  =  £(Pidt  +  Vi* x) 

i  =  l 


with  the  definitions  of  d^  and  given  in  Fig.  4a  as  the  projection  of  the 
geometric  path  onto  the  local  Cartesian  coordinates.  The  arrival  time  can  be 
determined  as  before  with 


defining 


Pi  =  Poj'Pz  =  Poj 


etc.  Thus, 


-  12 


sinQm  cos0„ 

Pm - - -  and  r)m  =  — — 

Om  am 


and  the  ~ — =  0  condition  leads  to  a  ray  going  from  the  source  to  the  receiver. 
The  first-motion  approximation  becomes 


Si  <x  Re(7-,i(p.,)  r„(p.,)  Si 


(24) 


Spreading  for  the  point  source  solution  becomes  slightly  more  complicated 
than  in  the  flat  case,  but  allowing 


V-= 

v  p 


i=iPt 


1/2 


(25) 


results  in  SP  defined  by 


v  r  r)  i  dp* 


(26) 


reducing  to 


=  (4,/ A)"*. 


at  the  direct  arrival  time.  The  details  of  this  result  have  been  given  previously 
by  Hong  and  Helmberger  (1975).  Thus,  the  point  source  solution  for  the 
geometry  given  in  Fig.  4b  becomes 


$  =  J_| 


1 

V? 


*dm 


(27) 


fcV- 
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where 


/(p)  -  Tlg(p  i)Tg3(pe) 


i! : ,  MM 
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Numerical  evaluation  of  (27)  yields  the  geometric  result  but.  also,  retains  longer 
period  information  since  (dp/dt)  can  be  evaluated  along  the  deHoop  contour  in 
the  usual  manner. 


Following  the  results  of  the  previous  section,  and  inserting  the  radiation 
pattern  for  the  simple  SH  motions  from  a  dislocation  source,  see  Helmberger 
and  Malone  (1975)  the  displacements  can  be  written 


v(r,z,Q,t)  =  2 


4rrp„  sinA 


d  !  a  / ,  \ 


^-JZ7(f)^£/lJ(0.\.i)F;(O 


where 


Vi)=S^(p)^2p)C,%Ti(p) 

<=1  Vfi  dt 


and  the  summation  of  n  rays  is  required.  The  various  symbols  are  defined  below: 


v(r, 2,0,0  =  displacement  on  free  surface 


M„  =  moment 


p0  =  density 

D(t)  =  dislocation  history 

Z?(t)  =  far-field  time  function 

/4 j(0.X, <5 )  =  cos20cosAsind  -  1/  2sin20sinAsin2d, 

i42(0.A,5)  =  -sin  0cos\cos<5  -  cosOsinAcos  2d, 

0  =  strike  from  the  end  of  the  fault  plane 
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\  =  rake  angle 
6  =  dip  angle 

r  =  distance  between  source  and  receiver 


p  =  ray  parameter 


1  2 
77  =  V~P\ 


A  =  epicentral  distance  in  radians 


=  correction  for  earth  flattening 


/9  =  shear  velocity 

and  where  the  vertical  radiation  patterns  are  given  by 


SH*  =  3- * 


c  n  + 1  z  >  h 

s  =  -i  S<A 


1  2 


The  correction  for  point  source  spreading  is  defined  by 


7  =  ~\/  —  V  — 
*  ^  PT  Pi  , 


This  solution  is  similar  to  the  flat  case  and  we  can,  obviously,  construct  the  div¬ 
ing  ray  response  for  a  smoothly  varying  structure  by  summing  the  primary  rays 
as  discussed  in  Fig.  2.  We  can  then  use  this  result  to  check  the  disk  ray  solution 
which  can  be  obtained  by  replacing  (30)  by 


=  SH,(p)Z£-{p)Z{6p/6t) 

Vfi 
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where  the  sum  is  taken  over  the  p  (t)  curve  as  described  by  Wiggins  (1976). 
For  a  simple  turning  ray  problem 


§IL\  =  _ l _ 

6t  jr-r(p)| 


(33) 


where  r(p)  is  distance  reached  by  a  ray  defined  by  p,  see  Fig.  5.  Substituting 
(33)  into  (32)  and  evaluating  (29)  yields  a  simple  step  response,  as  discussed  by 
Chapman  (1976).  Essentially,  (33)  has  a  simple  square-root  singularity  at  r  = 
r(p),  and  rays  that  hit  the  surface  near  the  receiver  dominate  the  behavior. 
Since  p  varies  along  the  path,  we  must  define  which  p  to  use  in  the  evaluation  of 
(32).  The  proper  choice  is  the  starting  p  at  the  source  as  outlined  in  the  previ¬ 
ous  section.  Note  that  for  the  case  of  an  up-going  direct  ray  the  two  methods 
can  be  interpreted  in  a  similar  manner.  Only  one  ray  is  involv.  d  in  both,  and 
applying  the  first-motion  approximation  of  (30)  yields  (32)  where  the  extra  (2)  is 
produced  by  the  double  valued  nature  of  expression  (33).  Thus,  the  application 
of  WKBJ  theory  to  the  locally  dipping  problem  appears  to  be  essentially  the  same 
as  for  the  uniform  layered  problem.  We  trace  the  ray  through  a  stack  of  layers 
down  to  the  turning  region,  turn  it  around  analytically,  and  follow  it  to  the  sur¬ 
face  obeying  Snell's  law.  The  special  treatment  at  the  turning  point  removes  the 
nonlinear  ray  parameter  effects  of  the  homogeneous  layered  parameterization 
as  is  well  known.  The  C,  factor  can  be  assumed  to  be  one  for  most  applications 
of  gentle  dipping  structure,  as  discussed  in  the  next  section  and  was  omitted 
from  (32). 

The  approach  followed  here  is  similar  to  that  followed  by  Wiggins  (1976)  in 
that  the  main  justification  for  expression  (33)  is  that  it  yields  results  compar¬ 
able  to  GRT.  A  theoretical  justification  of  applying  WKBJ  to  laterally  varying 
structure  is  given  by  Chapman  and  Drummond  (1982). 


3.  Applications 


In  this  section  we  will  briefly  outline  possible  applications  of  these  approxi¬ 
mate  solutions  to  seismological  problems.  First,  the  direct  or  up-going  energy 
problem  is  discussed  when  motions  in  the  sloping  tayers  of  a  sedimentary  basin 
are  excited  by  a  line-source.  In  this  form  Finite  Difference  calculations  can  be 
used  to  check  the  accuracy  of  the  GRT  results.  Next,  the  point  source  excitation 
of  Love  waves  is  considered  in  the  presence  of  sloping  structure  followed  by 
models  of  growing  Lithosphere.  Finally,  we  construct  synthetics  for  laterally 
varying  upper  mantle  models  and  confirm  the  usefulness  of  WKBJ  at  long 
periods. 


Local  seismograms 


One  of  many  complexities  associated  with  strong  motion  seismology  is  the 
noticeably  long  duration  of  high  frequency  P-waves  observed  in  sedimentary 
basins.  These  waves  are  generally  polarized  onto  the  vertical  component  due  to 
the  strong  velocity  gradients  near  the  surface.  The  latter  portion  of  these 
observed  motions  are  generally  depleted  at  lower  frequency.  Thus,  one  might 
conclude  that  there  are  propagational  waveguides  that  preferentially  prolong 
high  frequency  motions.  Non-planar  surface  layering  appears  to  have  this  pro¬ 
perty.  This  calculation  will  be  done  with  SH-waves  since  this  type  of  motion  is 
studied  throughout  the  remainder  of  this  paper,  but  we  would  expect  that  P- 
waves  would  behave  in  a  similar  manner. 

A  single  low  velocity  layer  which  grows  with  distance  between  the  source 
and  receiver  is  assumed  with  a  line  source  of  SH  motion  situated  at  a  depth  of 
5.5  km.  The  response  build-up  as  a  function  of  the  number  of  multiples  is 
displayed  in  Fig.  6.  The  square-root  singularity  indicated  in  expression  (13)  is 


apparent  for  the  direct  arrival.  Note  that  after  one  bounce  the  reflection  from 


the  lower  interface  becomes  complex  because  of  the  local  ray  parameter  effect 
and  a  head  wave  and  post-critical  angle  reflection  develops.  After  two  bounces, 
the  time  separation  between  the  head  wave  onset  and  reflection  times  becomes 
less  and  the  reflected  spike  increases  in  strength.  After  many  bounces  the  ray 
can  no  longer  reach  critical  angle  and  still  fit  into  the  waveguide.  Thus,  (R  )n 
becomes  small  since  the  reflection  coefficient  (R)  becomes  less  than  one.  The 
drop-off  in  amplitude  of  the  multiples  occurs  abruptly  at  this  time  on  the 
record. 


The  corresponding  point  source  response  displayed  in  Fig.  7  can  be 
obtained  from  expression  (29).  Neglecting  the  factor  produces  a  similar 
response  with  a  slight  reduction  in  later  arrivals,  roughly  13%  for  the  last  arrival. 
Thus,  point  source  amplitudes  can  be  approximated  quite  well  by  scaling  line 
source  results  by  the  square-root  of  the  distance  factor  similar  to  the  flat  case. 
Note  that  the  Cagniard-de  Hoop  technique  proves  particularly  useful  in  tracing 
these  rays  and  evaluating  their  individual  contributions.  However,  as  mentioned 
earlier,  this  series  of  rays  does  not  necessarily  converge  to  the  exact  solution 
and  some  demonstration  of  accuracy  is  required.  This  was  attempted  earlier  by 
Hong  and  Helmberger  (1976)  but  not  very  convincingly.  A  much  more  rigorous 
comparison  is  presently  being  conducted  by  Vidale  et  al  (1984)  with  respect  to 


using  line-source  numerics  to  construct  point-source  synthetics.  Preliminary 
results  of  the  comparison  of  the  two  techniques  for  this  single  model  is  displayed 
in  Fig.  8.  The  top  trace  is  the  broadband  result  displayed  in  Fig.  6.,  with  a 


filtered  response  in  the  middle  for  comparison  with  finite-difference  results  on 
the  bottom.  The  highest  frequencies  have  been  removed  in  this  comparison  due 
to  computational  expense  but  the  existence  of  strong  high  frequency  multiples 
is  striking.  Since  the  finite-difference  calculation  can  be  performed  on  any  arbi¬ 
trary  two-dimensional  structure  we  have  extended  the  thin  layer  directly  above 


the  source  to  the  left  as  a  flat  thin  layer  avoiding  the  wedge  effect  which  is  obvi¬ 
ously  omitted  in  the  ray  solution.  Comparison  with  and  without  the  wedge  and 
many  other  complexities  involving  double-couple  solutions  constructed  by  line- 
to-point  source  operators  are  discussed  in  Vidale  et  al  (1984).  We  will  suppose 
throughout  the  remainder  of  this  paper  that  the  generalized  ray  modifications 
discussed  in  the  previous  section  are  sufficiently  accurate  to  test  the  WKBJ  syn¬ 
thetics. 


Love  waves  at  Regional 


stances 


Another  interesting  application  of  the  above  technique  is  in  the  develop¬ 
ment  of  Love  waves  and  the  effects  of  traveling  across  oceanic-to-continental 
transitions.  This  problem  was  encountered  in  a  recent  paper  by  Grand  and 
Kelmberger  (1984b)  when  the  so-called  G-phase,  the  name  applied  to  the  impul¬ 
sive  Love  waves  associated  with  oceanic  paths,  interferes  with  mantle  arrivals. 
Apparently,  this  situation  occurs  for  well  developed  Lithospheres  associated  with 
older  plates  over-lying  slower  upper-mantle  velocity  models.  The  beginning  por¬ 
tions  of  the  G-phaso  as  recorded  slightly  inland  develop  longer  periods  than 
observed  at  Island  stations.  Their  period  and  arrival  times  are  compatible  with 
the  model  presented  in  Fig.  9a.  A  dipping  model  with  arrival  times  compatible 
with  the  flat  model  is  presented  in  Fig.  9b  along  with  th§  comparison  of  step 
responses  given  on  the  right.  Note  that  the  first  30  sec  of  motion  are  neariy 
identical.  The  higher  frequency  portions  of  the  Love  wave  become  less  pro¬ 
nounced  in  the  dipping  case  but  the  general  appearance  is  similar  to  the  pure- 
oceanic  case,  see  Grand  and  Kelmberger  (1984b). 


It  appears  that  as  the  Lithosphere  ages  it  gets  thicker,  for  example,  see 
Forsyth  (1975).  A  preliminary  model  of  predicted  Love  waves  for  this  situation  is 
given  in  Fig.  10,  also  included  are  synthetics  for  a  fast  and  slow  mantle.  The  long 
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period  nature  of  the  synthetics  from  the  dipping  model  is  similar  to  the  slow 
model  as  we  might  expect.  However,  there  is  considerable  roughness  at  the 
start  of  the  Love  waves  caused  by  the  mixed  paths  involving  both  the  crust  and 
lid. 

Observationally,  we  see  upper-mantle  arrivals  starting  near  these  ranges. 
Thus,  the  diving  energy  must  be  added  to  these  synthetics  following  the  strategy 
discussed  earlier.  This  can  be  accomplished  by  summing  GR’s  or  by  applying 
WKBJ. 


Upper  -mantle  models 

In  this  section  we  investigate  effects  of  lateral  variation  in  upper  mantle 
models,  as  displayed  in  Fig.  11.  We  have  chosen  a  particularly  simple  case  with 
no  low  velocity  zone  to  simplify  the  comparison  of  GRT  with  WKBJ  synthetics.  A 
further  simplification  is  made  by  allowing  the  two  models  to  be  connected  in  a 
linear  fashion  as  displayed  in  the  middle  column. 

Following  the  WKBJ  approach  we  first  illuminate  the  model  by  tracking  a  set 
of  rays  from  the  source  towards  prospective  receivers.  These  rays  reach  the 
surface  at  r(p)  in  time  T(p).  The  travel  time  at  a  particular  receiver,  r,  can  be 
written  t(p)  =  p  r  +  T(p)  -  p  r  (p).  Note  that  p  changes  in  each  layer  but  they  are 
all  functions  of  the  beginning  p.  Thus,  we  can  construct  the  t  versus  p  curves  as 
displayed  in  Fig.  12  for  reversed  profiles.  The  largest  ray  parameter,  pmitx,  is  .26 
which  corresponds  to  the  crustal  velocity  of  3.9.  Next,  we  perform  the  numeri¬ 
cal  derivative  ( 6p/6t )  of  these  curves.  Note  that  there  will  be  a  large  trunca¬ 


tion  phase  at  the  near  stations  at  pmtx.  This  can  be  avoided  by  including  the 
product  of  the  transmission  coefficients,  TC's,  across  the  Moho,  the  reference 


I  * 
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interface  discussed  earlier,  since  TC (pmax)  is  small.  Thus,  the  product  of  the 
TC's  with  ( 6p/6t )  has  a  relatively  smooth  behavior.  The  head  wave  along  the 
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Moho  is  added  in  by  including  the  reflected  generalized  ray.  By  performing  the 
convolution  indicated  in  expression  (29)  we  derive  step  responses  from  (t  vs.  p  ) 
curves  displayed  in  Fig.  12.  These  results  are  shown  in  Fig.  13.  Short  period 
synthetics  are  included  to  emphasize  the  rapid  decay  of  amplitude  at  the  tripli¬ 
cation  tips.  Eliminating  the  truncation  phase  discussed  here  can  also  be 
achieved  by  a  modification  of  the  Gaussian  beam  technique  as  developed  by 
Madariaga  and  Papadimitriou  (1984). 

The  synthetics  at  the  smallest  ranges  are  completely  controlled  by  the  shal¬ 
low  structure  and  the  local  model.  Thus,  the  first  arrival  from  the  Fast-to-Slow 
synthetics  have  a  shorter  travel  time  which  causes  the  triplication  from  the  40 0 
km  discontinuity  to  arrive  later  than  in  the  reverse  profile. 

A  more  detailed  plot  of  the  Slow-to-Fast  profile  is  displayed  in  Fig.  14  along 
with  the  GRT  responses  for  comparison.  The  synthetics  are  appropriate  for  the 
WWSSN  long  period  system.  A  typical  strike-slip  source  was  assumed  with  a  tri¬ 
angular  time  history  of  (1,  1,  1  secs)  and  a  t  *  =  3,  see  Grand  and  Helmberger 
(1984a). 

Note  that  there  is  a  distinct  change  in  the  latter  portion  of  the  WKBJ  step 
responses  between  17  and  18°.  This  is  caused  by  omitting  the  head  waves  from 
along  the  top  of  the  model  for  distances  beyond  17°.  However,  no  apparent 
change  in  the  synthetics  occurs  at  this  range  suggesting  that  the  long  period 
drift  is  outside  the  pass-band  of  the  operators  used  in  generating  these  synthet¬ 
ics.  The  high-frequency  spikes  so  apparent  in  the  GRT  step  responses  are  like¬ 
wise  removed  by  the  convolution  operators. 


The  400  km  discontinuity  is  treated  slightly  differently  in  the  two  methods 
which  leads  to  some  shifts  in  the  triplication  position.  In  GRT,  the  400  km 
discontinuity  is  generally  treated  as  a  sharp  jump  in  velocity  since  this  leads  to 
the  best  results  when  compared  against  reflectivity,  see  Burdick  and  Orcutt 
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(1978).  On  the  other  hand.  WKBJ  requires  a  smooth  transition,  3  km  transition  in 
this  particular  case,  such  that  the  (p  vs.  t  )  curve  is  smooth.  Thus,  the  sharp 
spikes  occurring  in  the  GRT  step  responses  near  14°  are  precritical  angle 
reflections  from  the  400  km  discontinuity.  Similarly,  the  triplication  seems  to 
extend  to  greater  distances  in  the  GRT  results.  Note  that  the  most  severe 
mismatch  occurs  near  this  range.  At  larger  ranges  the  two  methods  agree  quite 
well,  especially  the  synthetic  waveforms.  In  fact,  the  synthetic  waveforms  agree 
at  all  distances  with  the  maximum  deviation  in  amplitude  of  about  25%.  And, 
since  these  synthetic  waveforms  are  used  to  interpret  observations  which  can 
seldom  be  modeled  as  well  as  the  agreement  between  these  two  methods,  we  can 
consider  the  WKB  modifications  successful.  For  more  precision  involving  sharp 
boundaries  we  suggest  breaking  the  p  integration  into  a  combination  of  WKB  for 
the  smooth  portion  of  the  model  and  a  generalized  ray  for  the  reflecting  inter¬ 
face,  for  example  see  Given  (1984). 

Conclusions 

In  this  paper  we  presented  a  hybrid  procedure  of  generating  complete 
seismograms  in  laterally  varying  structure  by  applying  a  mixture  of  GRT  and 
DRT.  First,  we  reviewed  the  modifications  of  GRT  required  for  dipping  structure 
in  terms  of  local  coordinates  and  ray  parameter  concepts  for  line  and  point 
source  theory.  Solutions  calculated  by  this  approach  not  only  agree  with 
geometric  results,  but  also  agree  with  longer  period  motions  such  as  computed 
with  finite-difference  methods.  Using  the  correspondence  between  GRT  and  DRT, 
disk  ray  theory,  we  can  express  the  latter  in  relatively  simple  form,  essentially 
applying  a  square-root  of  distance  correction  to  line  source  spreading.  Com¬ 
parisons  between  GRT  and  DRT  synthetics  of  diving  energy  paths  agree  reason¬ 
ably  well.  Thus,  we  can  construct  nearly  complete  seismograms  with  a 
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combination  of  GRT  and  DRT  with  the  former  used  to  handle  the  shallow  struc¬ 
ture.  Some  useful  demonstrations  of  the  methods  are  given  for  crustal  and 
upper-mantle  models. 
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Abstract 

Ve  have  modeled  the  SH  action  from  earthquakes  in  the 
northwest  Atlantic  ocean  to  derive  an  oceanic  upper  mantle  shear 
velocity  model.  The  signals  were  recorded  on  long-period  VUSSN  and 
Canadian  network  stations  on  the  east  coast  of  North  America.  The 
travel  times  and  waveforms  of  seismograms,  in  the  distance  range  11°  to 
16°,  were  used  to  constrain  the  lid  structure.  This  data  Indicates  a 
fast  (4.75  km/sec)  lid  of  about  100  km  thickness  in  the  older  western 
Atlantic.  From  16°  to  28°  the  S  waves  pass  through  triplications  due  to 
velocity  jumps  near  400  km  and  660  lun.  The  branches  from  both 
discontinuities  are  visible  in  the  S  wave  data.  These  arrivals  were 
modeled  using  synthetic  seismograms  to  obtain  accurate  travel  times. 
The  two  triplications  are  also  apparent  in  SS  data  from  30°  to  50°.  The 
triplication  curves  derived  from  the  S  waves  agree  with  the  SS  data. 
Given  the  lid  structure,  the  waveforms  and  traveltimes  from  the  more 
distant  data  put  tight  constraints  on  the  shear  velocities  at  greater 
depths.  The  velocity  below  200  km  was  found  to  be  indistinguishable 
from  a  model  of  the  East  Pacific  Rise  (Grand  and  Helmberger,  1983)  found 
using  the  same  technique.  This  is  in  contrast  to  the  Canadian  shield 
model  found  in  the  same  paper.  We  find  the  Canadian  shield  to  be  faster 
than  both  the  old  northwestern  Atlantic  and  the  young  East  Pacific  Rise 
to  about  400  km  depth.  No  variations  below  400  km  are  necessary  to 
explain  the  data. 


Introduction 


Upper  mantle  shear  velocity  has  been  shown  in  numerous  studies  to 
vary  significantly  under  different  tectonic  provinces.  The  sensitivity 
of  shear  velocity  to  temperature  makes  the  study  of  lateral  variation  in 
shear  structure  an  important  tool  in  determining  the  thermal  state  of 
the  earth.  There  has  been  debate,  however,  on  the  depth  extent  of  shear 
heterogeneity  between  different  tectonic  regimes.  Jordan  (1981)  has 
proposed  that  old  cratons  have  deep  roots,  from  400  km  to  650  km  in 
depth.  This  idea  is  based  largely  on  the  differences  in  time  of  ScS 
multiples  in  different  regions  (Sipkin  and  Jordan  1976,1980).  Anderson 
(1979),  on  the  contrary,  states  that  continental  roots  extend  no  deeper 
than  about  150  km  to  200  km.  Okal  and  Anderson  (1975)  and  Anderson 
(1979)  claim  that  the  ScS  data  and  surface  wave  data  are  incompatible 
with  the  'deep  root'  hypothesis  of  Jordan.  It  will  take  detailed 
studies  of  the  upper  mantle  beneath  several  cratons  and  many  oceanic 
areas  before  generalizations  about  the  deep  structure  of  continents  and 
oceans  can  be  stated  with  certainty. 


Grand  and  Helmberger  (1983),  paper  1,  studied  the  upper  mantle 
under  North  America  and  part  of  the  East  Pacific  Rise  using  S  and  SS 
bodywave  phases  in  the  distance  range  10°  to  60°.  They  showed  that 
using  the  travel-times  and  waveforms  of  S  and  SS  one  could  constrain  the 
shear  structure  to  800  km  depth  with  a  resolution  of  about  1%  over 
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100  km.  This  technique  provides  better  resolution,  at  depth,  than 
conventional  surface  wave  dispersion  studies  due  to  the  fact  that  the 
wavelengths  of  energy  sampling  deeper  structure  do  not  increase 
substantially.  In  this  paper  we  have  studied  the  western  part  of  the 
North  Atlantic  ocean  using  essentially  the  same  technique  of  modeling  S 
and  SS  phases.  The  area  studied  has  oceanic  crust  ranging  in  age  from  0 
to  150  m.y.  (Sclater  &  Parsons,  1981).  Thus  we  have  an  opportunity  to 
compare  the  upper  mantle  structure  of  an  active  ridge,  old  ocean, 
tectonic  continent  and  stable  craton  derived  from  a  single  method. 

The  data  from  the  Atlantic  ocean  show  features  very  similar  to  the 
data  from  the  East  Pacific  Rise  and  the  Canadian  shield. 
Discontinuities  near  400  km  and  660  km  produce  triplication  branches 
which  are  visible  in  both  the  S  and  SS  waves.  Using  the  times  from 
these  branches  we  have  derived  an  average  model  for  the  older  part  of 
the  Atlantic  (70  to  150  m.y.)  near  North  America.  The  model,  listed  in 
table  1,  and  shown  in  figure  1  with  the  previously  derived  models  from 
paper  1,  has  a  high  velocity  lid  about  100  km  thick.  Below  the  lid,  the 
velocity  decreases  to  200  km  depth  at  which  depth  it  has  the  same 
velocity  as  that  derived  for  the  East  Pacific  Rise  in  paper  1.  Below 
200  km  we  find  no  difference  in  structure  between  the  old  Atlantic  and 
the  young  East  Pacific  regions. 
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To  determine  the  upper  mantle  shear  velocity  beneath  the  Atlantic, 
we  modeled  the  SH  motion  from  Atlantic  earthquakes  recorded  on  long 
period  WWSSN  (World  Wide  Standard  Seismic  Network)  and  CSN  (Canadian 
Seismic  Network)  stations  in  eastern  North  America.  The  approach  is  to 
match  the  travel  times  and  waveforms  of  nearin  data,  using  synthetic 
seismograms,  to  constrain  the  shallow  structure.  As  one  looks  at 
seismograms  at  further  distances,  deeper  structure  is  sampled.  Using 
synthetic  seismograms  we  can  analyze  what  from  depth  individual  arrivals 
come  and  thus,  adjust  the  model  appropriately  when  discrepancies  between 
the  synthetics  and  data  appear.  The  starting  model  was  taken  as  the 
East  Pacific  Rise  model  (TNA)  found  in  paper  1.  Changes  to  it  were  made 
from  the  top  downward  until  the  Atlantic  data  was  satisfied.  The 
resolution  of  this  trial  and  error  approach  was  investigated  by  Given 
(1983).  He  developed  a  formal  inversion  for  structure,  using  S  wave 
travel  times  and  waveforms,  and  applied  it  to  data  from  the  East  Pacific 
Rise  and  western  North  America.  The  resulting  model  was  very  similar  to 
TNA,  which  was  derived  by  a  trial  and  error  technique  thus  we  feel 
confident  in  the  conclusions  of  our  modeling. 

The  earthquakes  used  are  listed  in  table  2.  Figure  2  shows  the 
locations  of  the  events  and  stations  relative  to  oceanic  crustal  age 
(Jordan  1981).  The  events  were  selected  on  the  basis  of  the  quality  of 
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the  SH  signals  from  them.  The  stations  used  were  nearly  naturally 
rotated  and  thus  data  with  strong  SV  energy  was  apparent.  Due  to  the 


problem  of  shear-coupled  PL  interference,  discussed  in 
Helmberger  &  Engen  (1974),  we  rejected  any  station  which  had  a  high  SV 
to  SH  ratio  at  the  arrival  time  of  S  or  SS.  The  events  on  the 
Mid-Atlantic  Ridge  had  teleseismic  SH  waves  similar  to  those  from 
stike-slip  earthquakes,  thus  we  assumed  the  mechanism  to  be  strike-slip 
for  them.  Event  1,  near  Bermuda,  was  studied  by  Stewart  and  Helmberger 
(1981),  and  event  2,  by  Molnar  and  Sykes  (1967).  The  focal  mechanisms 
of  the  two  events  were  thrust  but  they  still  produced  simple  SH  waves  to 
some  stations. 

The  data  fall  in  the  range  11°  to  52°  and  sample  ocean  of  varying 
age.  From  11°  to  30°,  the  S  waves  travel  wholly  in  the  upper  mantle  and 
the  triplications  are  obvious  in  the  data.  Data  at  these  distances  came 
from  events  1  and  2,  and,  as  can  be  seen  on  the  map,  their  paths  were 
totally  in  ocean  older  than  100  m.  y..  Data  from  27°  to  52c  came  from 
the  events  near  the  Mid-Atlantic  Ridge.  We  used  the  SS  waveforms  and 
SS-S  times  at  these  ranges  to  determine  the  structure.  SS  goes  through 
the  same  triplication  as  S  at  nearer  distances  as  discussed  in  paper  1. 
These  earthquakes  are  in  much  younger  ocean  than  those  used  for  S, 
however,  their  midpoints  are  all  within  ocean  older  than  100  m.  y.. 
Using  SS-S  times  should  decrease  the  contamination  of  the  data  by  the 
initial  propagation  through  younger  structure. 

To  derive  an  upper  mantle  structure  which  would  reproduce  the 
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relatively  complicated  seismograms  used,  we  constructed  synthetic 
seismograms.  Synthetics  enable  one  to  model  interfering  signals  and, 
also,  to  intuitively  understand  them.  The  procedure  used  to  calculate 
the  synthetics  is  well  established.  One  assumes  a  seismogram  consists 
of  a  series  of  linear  operators  representing  the  different  factors  in  a 
record  (Helmberger  &  Burdick,  1978).  We  can  write  a  seismogram  Y(t)  as 

Y(t)-S(t)*I(t)*A(t)*M(t) 

where  S(t)  is  the  source,  l(t)  is  the  instrument  response,  A(t)  is 
an  attenuation  operator  and  M(t)  is  the  Green's  function  for  wave 
propagation  through  the  mantle. 

The  source  function  was  constructed  by  fitting  the  teleseismic  S 
waveforms.  A  trapezoidal  time  function  was  convolved  with  an  operator 
representing  the  S,sS  interaction,  dependent  on  focal  mechanism 
(Langston  &  Helmberger,  1975).  The  depth  and  time  function  are  adjusted 
to  fit  the  teleseismic  data.  All  the  events  used  in  this  study  are 
shallow  and  produced  simple  teleseismic  SH  waves.  The  details  of  the 
source  are  not  too  important  as  the  arrivals  we  modeled  are  generally 
well  separated. 


The  attenuation  operator  was  taken  from  Futterman  (1962)  and 
Carpenter  (1965).  The  amount  of  attenuation  is  parameterized  by  the 
value  of  t*,  the  travel-time  divided  by  the  average  Q  along  the  ray 
path.  It  is  generally  believed  that  Q  is  lower  in  the  upper  few  hundred 
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km  than  below.  For  this  reason,  we  attenuated  SS  arrivals  more  than  S 

it  ^ 

arrivals.  At  of  4  was  used  for  S  and  at  of  7  for  SS  waves.  These 
values  are  rather  arbitrary  but  are  fairly  close  to  values  computed  from 
published  Q  models  such  as  PREM  (Dzeiwonski  &  Anderson,  1981). 
Improvements  in  Q  structure  will  help  the  synthetic  to  data  amplitude 
fits  but  should  not  affect  the  relative  timings  of  different  arrivals 
significantly. 


The  Green's  functions,  for  upper  mantle  propagation,  were  computed 
using  two  techniques.  Preliminary  modeling  of  the  data  was  done  using 
the  WKBJ  method  of  Chapman  (1978)  and  Wiggins  (1976).  This  method  is 
very  inexpensive  and  useful  in  the  modeling  process.  There  are 
situations,  however,  where  it  does  a  poor  job,  such  as  energy  from 
shadow  zones  and  interactions  with  sharp  discontinuities.  Therefore, 
after  deriving  a  model  fitting  the  traveltimes  of  the  most  obvious 
arrivals  with  the  WKBJ  synthetics,  we  used  a  Cagniard-de  Hoop  code 
(Helmberger,  1973)  to  model,  in  more  detail,  our  data.  Comparisons  of 
synthetics  by  the  two  methods  are  given  in  paper  1  for  two  upper  mantle 
models . 


The  Cagniard-de  Hoop  synthetic  is  generated  by  adding  the  responses 
of  generalized  rays  through  a  stack  of  homogeneous  layers.  A  limitation 
of  the  technique  is  that  for  SS,  or  other  multiple  bounce  phases,  a 
large  number  of  rays  are  needed  to  compute  the  response.  At  different 
distances  we  empirically  found  which  rays  were  necessary  by 
progressively  adding  more  rays  until  the  synthetics  did  not  change 


appreciably.  The  rays  used  in  different  situations  will  be  mentioned  in 
the  following.  Comparisons  of  Cagniard-de  Hoop  with  the  reflectivity 
method  (Burdick  &  Orcutt,  1978)  have  been  positive  and  thus  we  feel  this 
approach  is  valid. 


Lid  Structure 


To  investigate  the  upper  mantle  below  200  km  a  good  estimate  of  the 
upper  200  km  is  necessary  since  all  the  data  used  to  look  at  deeper 
structure  will  be  contaminated  by  the  more  shallow  structure.  The  upper 
part  of  the  model  was  found  by  fitting  the  traveltimes  and  waveforms  of 
the  SH  motion  at  closein  distances  (11°  to  16°).  For  a  continental 
structure,  at  these  ranges,  there  are  well  separated  arrivals  from  the 
mantle.  Depending  on  the  structure,  they  can  be  from  above  or  below  the 
lid.  These  separate  arrivals  (S,SS,SSS  etc.)  were  modeled  in  paper  1 
to  determine  the  average  shallow  structure  beneath  the  Canadian  shield 
and  western  North  America.  The  reason  the  arrivals  are  well  separated 
at  these  ranges  is  due  to  the  thick  crust  above  the  mantle.  Each  time  a 
phase  comes  to  the  surface  it  must  travel  through  the  much  slower  crust 
and,  thus,  multiple  bounce  phases,  such  as  SS  and  SSS,  fall  rapidly 
behind  the  direct  arrival  S  wave.  In  an  oceanic  structure  the  crust  is 
much  thinner  and  thus  the  separation  of  multiple  bounce  arrivals  is  much 


less  in  time.  In  figure  3  we  illustrate  this  effect  using  the 
Cagniard-de  Hoop  generalized  ray  technique.  We  show  the  effect  of 
adding  successive  multi-bounce  arrivals  for  two  simple  structures.  To 
represent  an  oceanic  structure,  we  used  a  7  km  thick  layer  of  velocity 
3.7  km/sec  over  a  half-space  with  velocity  4.75  km/sec.  The  continental 
synthetics  are  for  a  30  km  layer,  of  velocity  3.7  km/sec  over  the  same 
half-space.  The  rays  used  are  shown  to  the  left  of  the  synthetic  step 
responses.  The  seismograms  are  calculated  at  a  distance  of  8°.  The 
front  of  the  records  are  due  to  multiple  bounce  head  waves  while  the 
high  frequency  later  arrivals  are  reflections.  Each  multiple  bounce  ray 
is  very  similar  in  the  two  structures,  however,  the  difference  in 
crustal  thickness  changes  the  time  separation  of  the  arrivals  and,  thus, 
the  interference.  This  produces  the  vastly  different  dispersion  in  the 
two  structures.  The  final  panel  is  the  response  after  convolving  with  a 
long-period  WWSSN  instrument  and  a  6  sec  trapezoidal  time  function. 
Obviously,  the  oceanic  data,  at  nearin  distances,  will  have  a  very 
different  appearance  than  the  continental  data  just  due  to  the  crust. 
Large  sets  of  multiply  reflected  crustal  rays  must  be  used  to  model 
oceanic  records  at  upper  mantle  distances.  In  terms  of  modes,  this  can 
be  understood  by  the  fact  that  the  fundamental  Love  mode,  from  20  to  60 
sec  travels  at  mantle  velocities  and  not  crustal  velocities,  as  in  the 
continent.  This  phenomena  was  pointed  out  by  Thatcher  and  Brune  (1969) 
in  their  discussion  of  the  interference  of  higher  mode  Love  waves  with 
the  fundamental  at  periods  less  than  60  sec  for  oceans. 


Inspite  of  the  fact  that  we  cannot  model  nearin  data  as  distinct 
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mantle  arrivals,  the  waveforms  still  contain  much  information.  Figure  4 
illustrates  how  we  might  use  data  at  nearin  distances,  say  12°. 
Figure  4a  shows  the  effect  of  changing  crustal  thickness.  The 
synthetics  were  computed  summing  rays  with  up  to  40  multiples  in  the 
crustal  layer.  A  strike-slip  source  is  at  5  km  depth.  With  changing 
crustal  thickness  the  arrival  time  of  the  first  downswing  stays  the 
same.  This  is  basically  headwave  energy  with  only  a  few  free  surface 
reflections.  The  later  swings,  in  the  waveform,  change  dramatically 
with  crustal  thickness  differences  of  only  a  few  kilometers.  The  signal 
is  higher  frequency  with  thinner  crust  due  to  the  fact  that  the  later 
arrivals  in  the  wavetrain  have  bounced  many  times  in  the  crust. 
Changing  the  crustal  thickness  changes  the  time  separation  of  the 
multiply  reflected  arrivals  and  this  becomes  more  important  in  the  back 
of  the  waveform  where  the  reflections  are  more  numerous.  Changing 
crustal  velocity  will  have  the  same  effect.  A  slower  crust  would  be 
equivalent  to  a  thicker  crust. 

Figure  4b  shows  the  effect  of  changing  lid  velocity  but  keeping  the 
crust  constant.  The  waveforms  are  not  changed  significantly  but  the 
absolute  travel-time  changes  by  12  sec  with  a  change  in  lid  velocity 
from  4.6  km/sec  to  4.8  km/sec.  In  figure  4c  we  have  added  a  further 
variation  by  putting  a  lower  velocity  half-space  below  the  lid.  In  this 
case,  all  the  crustal  reverberations  were  used  along  with  rays  which 
have  reflected  from  the  bottom  of  the  lid  up  to  4  times.  The  reflection 
coefficient  at  the  lid  half-space  boundary  is  the  opposite  sign  of  that 
at  the  crust  lid  boundary.  Thus,  for  a  thin  lid,  the  lid  rays  cancel 


the  front  of  the  waveform.  The  thicker  lids  have  little  effect  at  12° 
because  the  rays  arriving  at  the  station  from  deeper  depths  are  steeper 
and  thus  the  reflection  coefficients  are  fairly  small.  The  thicker  lids 
will  affect  the  Love  wave  at  further  distances. 

We  approached  the  data  at  nearin  distances  in  the  following  way. 
First,  the  crustal  thickness  was  adjusted  to  match  in  time  the 
successive  peaks  and  troughs  in  the  wavetrain.  The  crustal  velocity  was 
kept  fixed  at  3.7  km/sec.  Next  the  lid  velocity  was  modified  to  match 
the  absolute  travel-time  of  the  first  downswing.  There  is  some  tradeoff 
here  between  the  lid  thickness  and  the  lid  velocity.  However,  as  can  be 
seen  in  figure  4,  there  is  a  difference  between  a  high  velocity  thin  lid 
and  a  low  velocity  thicker  lid.  The  separation  of  the  first  two 
downswings  is  different  in  the  two  cases  since  lid  velocity  basically 
shifts  the  whole  waveform  but  lid  thickness  controls  only  the  front  of 
the  signal.  The  lid  models  found  in  this  way  are  obviously  very  crude 
but  they  should  constrain  the  main  features  of  the  very  shallow  mantle 
and  allow  us  to  investigate  deeper  structure. 

Figure  5  shows  the  two  nearest  seismograms  used  in  the  study.  The 
later  arrivals  at  BEC  were  too  large  to  recover.  These  seismograms  both 
recorded  almost  pure  SH  motion  with  paths  totally  in  ocean  older  than 
100  m.  y..  The  dispersion  of  the  signal  at  SJG  required  a  crustal 
thickness  of  9.5  km.  This  is  thicker  than  normal  oceanic  crust,  though 
it  could  be  that  we  used  too  high  a  crustal  velocity.  A  low  velocity, 
low  Q  layer  at  the  surface  would  decrease  the  amplitude  of  the  back  of 
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the  wavegroup,  since  the  signal  there  is  due  to  more  surface  reflections 
than  the  front.  It  would  also  decrease  the  average  crustal  velocity. 
Officer  et  al.  (1952)  report  an  average  2  km  thick  sediment  layer  south 
of  Bernuda  from  a  refraction  survey.  Another  cause  of  this  could  be 
that  we  used  the  wrong  depth  for  the  source.  The  source  depth  for  the 
synthetics  is  5  km  but  the  event  could  be  deeper.  For  the  long  periods 
this  will  have  no  effect  but  the  shorter  periods  could  be  changed 
slightly.  To  match  the  timing  of  the  first  arrival  a  lid  velocity  of 
4.75  km/sec  was  needed.  The  separation  of  the  first  downswings  requires 
a  lid  thickness  of  at  least  80  km.  At  these  distances  a  thicker  lid 
would  appear  about  the  same.  Primary  rays  from  below  the  lid  are  not 
noticeable. 

In  figure  6  we  show  further  data  at  regional  distances.  These 
stations  are  within  continents  and  the  transition  from  oceanic  to 
continental  crust  has  an  obvious  affect  on  the  waveforms.  The  overall 
period  of  these  signals  is  much  longer  than  the  SJG  record  from  the  same 
event.  From  the  previous  discussion,  this  is  what  one  would  expect  for 
an  increase  in  crustal  thickness.  We  used  a  crust  of  19.5  km,  with  the 
same  velocity  as  before,  to  model  these  records.  This  was  found 
empirically  to  fit  the  overall  dispersion  of  the  data.  How  the  actual 
ocean-continent  transition  affects  the  Love  wave  is  not  clear  but  since 
the  crustal  thickness  does  not  change  the  arrival  time  by  much,  using  an 
average  thickness  for  the  crust  is  probably  sufficient.  Note  that  the 
arrival  times  are  well  fit  for  this  data  set  using  the  same  model  as  in 
the  previous  figure.  Since  we  have  fit  data  in  three  directions,  from 


the  same  event,  with  the  same  model,  it  appears  the  event  was  not 
greatly  mislocated,  nor  is  there  a  large  directional  anisotropy.  The 
arrows,  in  the  figure,  indicate  the  reflected  arrival  from  the  400  km 
discontinuity.  It  can  also  be  seen  in  the  data  at  about  the  right  time. 
The  data,  in  the  previous  figures,  constrain  the  lid  to  have  a  high 
velocity  (4.75  km/sec)  and  a  thickness  of  at  least  80  km.  The  timing  of 
the  400  km  reflection  indicates  the  velocities  below  the  lid  must  be 
significantly  lower  than  those  found  in  the  shield  in  paper  1.  Further 
evidence  for  the  thick  lid  will  be  presented  later  in  relation  to  the  SS 
data. 


Structure  below  the  lid 


Figure  7  shows  the  travel-time  curve  for  the  Atlantic  model.  The 
previous  data  basically  established  the  AB  branch  to  16°,  this  is  energy 
from  the  lid.  The  two  seismograms  at  16.3  and  16.8°  showed  an  arrival 
corresponding  to  the  cusp  labeled  C,  a  reflection  from  the  405  km 
discontinuity.  The  dashed  lines  refer  to  diffracted  energy  or  tunnelled 
energy  through  the  lid.  There  are  several  more  seismograms  from  events 
1  and  2  which  sample  deeper  than  the  lid.  Figure  8  shows  a  small 
profile  to  the  north  of  these  events.  Again,  most  of  the  paths  are 
through  old  ocean  but  the  receivers  are  on  continental  crust.  In  the 


synthetics  we  again  used  an  'average'  crustal  thickness  of  19.5  km.  The 
crustal  arrivals,  or  Love  waves,  arrive  just  after  the  branch  labeled  B. 
Errors  in  our  crustal  approximation  should  not  affect  the  measurements 
of  the  mantle  arrivals  which  determine  the  triplication  curve.  The 
mechanism  appropriate  for  event  2  was  used  in  the  synthetics.  This  is 
inappropriate  for  the  record  at  SCH  from  event  1  causing  some  of  the 
misfit  there,  though  there  are  obvious  arrivals  in  the  record 
corresponding  to  those  in  the  synthetic  at  the  right  time.  With  respect 
to  figure  7,  we  can  follow  the  three  branches  over  several  degrees. 
Branch  F  is  from  below  the  660  km  discontinuity,  branch  D  is  a  wide 
angle  reflection  from  the  660  and  branch  B  is  the  wide  angle  reflection 
from  the  405  km  discontinuity.  The  model  derived  from  East  Pacific  Rise 
data  in  paper  1,  which  we  used  below  200  km,  fits  the  timing  of  all  the 
branches  very  well.  Differences  of  1%  over  200  km  at  any  depth  range 
would  be  noticeable  in  the  relative  timings  of  these  phases. 

In  figure  9,  two  records  from  this  profile  are  shown  compared  to 
synthetics  computed  from  the  Atlantic  model  and  a  modified  Atlantic 
model.  The  modified  model  is  the  same  as  the  Atlantic  model  except  from 
225  to  400  km  depth,  where  we  have  substituted  the  velocities  found 
beneath  the  Canadian  shield  at  those  depths.  The  SCH  synthetics  were 
done  using  the  appropriate  source  for  event  1,  also  the  Love  wave  was 
not  Included  as  this  station  is  much  further  inland  than  the  others  and 
the  Love  wave  arrives  later.  This  figure  illustrates  the  sensitivity  of 
the  arrival  time  of  branch  B  to  the  structure  from  200  to  400  km.  The 
arrows  indicate  the  arrival  from  the  400  km  discontinuity  in  the  data 
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and  synthetics.  The  modified  model  is  too  fast  by  about  5  sec.  The 
velocity  could  be  decreased  from  100  to  200  km  to  compensate  for  this 
but  the  observations  of  branch  C  at  16°  and  in  the  SS  data  (shown  later) 
limit  the  amount  one  can  modify  the  structure.  Note  that  the  660 
triplication  is  hardly  changed  at  all  by  the  modification. 

From  this  small  data  set  we  can  conclude  several  things.  With  the 
constraint  of  a  high  velocity  lid,  the  velocity  from  200  to  400  km  must 
be  substantially  the  same  as  beneath  the  East  Pacific  Rise.  Also,  the 
structure  below  400  km  appears  to  be  identical  to  both  the  Canadian 
shield  and  the  East  Pacific  Rise. 


Observations  of  SS 


Events  3  through  9  are  further  from  the  stations  than  the  previous 
data,  but  they  can  still  be  used  for  studying  the  upper  mantle.  Events 
3,4  and  5,  in  particular,  have  paths  to  east  coast  stations  through 
ocean  mostly  older  than  100  m.  y. .  These  events  produced  data  in  the 
range  27  to  53°,  thus  the  SS  phases  from  these  earthquakes  are  going 
through  the  same  triplications  as  S  waves  from  13.5  to  26.5°.  As 
discussed  in  paper  1,  SS  data  in  this  range  can  be  very  useful  in 
modeling  the  mantle.  The  S  wave  bottoms  below  800  km  where  the 


structure  is  fairly  well  known.  The  source  and  receiver  structure  is 
the  same  for  both  S  and  SS,  thus  the  S  wave  serves  as  a  useful  absolute 
time  marker  for  the  SS  phase  which  travels  through  the  heterogeneous 
upper  mantle. 

Figure  10  illustrates  the  synthetic  construction  of  seismograms  in 
this  range.  These  seismograms  were  computed  using  the  Cagniard-de  Hoop 
technique.  The  first  column  shows  the  results  of  just  adding  in  primary 
S  rays.  The  triplication  branches  are  labeled  on  the  diagram.  The 
backbranch  due  to  the  660  discontinuity  can  be  followed  to  near  32°. 
The  400  backbranch  is  seen  to  about  30°  but  there  is  a  small  arrival 
continuing  to  greater  distances.  This  more  distant  phase  is  not 
actually  from  the  400  km  discontinuity  but  is  due  to  rays  coming  just 
out  of  the  shadow  zone,  bottoming  near  250  km  depth,  which  travel  to 
large  distances  because  of  the  high  velocity  thick  lid.  If  there  is  a 
low  Q  layer  below  the  lid  this  arrival  could  be  more  attenuated.  In  the 
second  column  we  have  added  what  we  call  the  Love  rays.  These  are 
multiples  in  the  lid  and  crust.  For  this  example,  a  19.5  km  crust  was 
used,  as  the  following  data  are  all  from  continental  stations.  Rays 
with  25  reflections  in  the  crust  and  up  to  4  reflections  from  the  bottom 
of  the  lid  were  used.  One  multiple  within  the  lid  was  also  found  to  be 
necessary.  One  can  see  that  this  energy  drastically  changes  the 
appearence  of  the  seismograms  near  the  400  km  backbranch  time.  As 
stated  before,  this  arrival  does  not  Interfere  with  upper  mantle 
arrivals  in  continental  areas.  Finally,  in  the  third  column,  we  have 
added  SS  rays  from  below  the  lid.  The  major  effect  of  these  rays  is 


seen  from  37  to  40°.  The  SS  B  branch  from  the  400  km  discontinuity  can 
be  seen  clearly.  Even  at  40°  it  is  well  behind  the  Love  wave. 


The  Love  wave,  at  these  distances,  again  provides  information  on 
the  thickness  of  the  lid.  In  figure  11  we  have  computed  just  the  Love 
rays  at  a  distance  of  40°.  The  crust  is  fixed  at  19.5  km  and  the  lid 
velocity  is  4.75  km/sec,  as  was  found  before.  We  have  varied  the  lid 
thickness  from  50  to  110  km.  Again,  the  rays  reflecting  from  the  bottom 
of  the  lid  are  negative  relative  to  those  of  the  crustal  rays.  The 
thinner  lids  have  stronger  reflections  from  this  boundary  and  thus  there 
is  more  destructive  interference  of  the  front  of  the  Love  wave.  Figure 
11  shows  the  response  of  a  WWSS  instrument,  convolved  with  a  5  sec 
trapezoidal  time  function,  to  the  Love  rays.  The  arrival  time  of  the 
first  large  upswing  is  very  sensitive  to  the  lid  thickness.  The  front 
of  the  Love  wave  from  27  to  40°  should,  then,  put  constraints  on  the 
lid.  The  other  arrivals,  at  these  distances,  will  constrain  the  deeper 
structure. 

Figure  12  shows  the  data  in  the  range  discussed  above  with  the  same 
synthetics  as  were  derived  in  figure  10.  The  traveltime  lines  in 
figure  10  are  again  reproduced  in  this  figure  on  the  data  and 
synthetics.  The  stations,  listed  next  to  the  seismograms  recorded  by 
them,  are  of  varying  distances  from  the  ocean.  We  used  a  19.5  km  thick 
crust  as  an  approximation  to  the  actually  more  complicated  situation. 
This  thickness  seems  to  satisfy  the  dispersion  of  the  later  part  of  the 
Love  waves  recorded  at  these  stations  and  should  not  affect  the  mantle 
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arrivals  too  much.  All  the  arrivals  pointed  out  in  the  previous 
discussion  of  the  synthetics  can  be  seen  in  the  data.  The  arrival  time 
and  amplitude  decay  of  the  660  km  discontinuity,  branch  D,  is  modeled 
extremely  well.  The  structure  near  600  km  was  derived  from  a  totally 
different  area  so  that  this  fit  is  quite  surprising.  This  indicates 
that  heterogeneity  below  400  km  from  the  East  Pacific  Rise  to  the 
Canadian  shield  to  the  old  northwest  Atlantic  is  extremely  small. 
Velocity  differences  of  1%  from  400  to  600  km  would  change  the 
synthetics  significantly  at  these  ranges. 

From  38  to  40°  the  B  branch  of  the  400  km  discontinuity  is  visible 
in  the  data,  as  an  SS  arrival.  It  can  be  identified  by  the  higher 
frequency  content  than  the  preceding  Love  wave  pulse  and  its  moveout 
from  the  S  wave.  Data  beyond  40°  confirm  our  identification  of  the  SS  B 
branch.  Note,  also,  that  the  timing  is  consistent  with  the  arrival  at 
16°,  figure  6,  which  we  modeled  as  the  400  km  reflection.  The  midpoints 
of  these  paths  are  in  ocean  older  than  100  m.  y.. 

The  part  of  this  section  showing  the  most  variability  in  the  data 
is  indicated  by  the  line  labeled  B.  The  arrivals  just  following  this 
line  are  due  to  rays  from  several  depths  as  was  illustrated  in 
figure  10.  The  front  part  of  this  waveform  is  due  to  S  rays,  and 
crustal  multiples  of  them,  coming  from  below  the  lid.  This  arrival  is 
most  sensitive  to  the  velocity  from  150  to  350  km.  The  late  time  of 
this  phase,  in  the  data,  indicate  low  velocities  at  this  depth  range. 
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Just  following  this  energy  is  the  start  of  the  Love  wave  which  travels 


at  near  the  lid  velocity.  The  time  of  this  wave  group  indicates  a  thick 
fast  lid  as  discussed  earlier.  This  energy,  however,  varies  by  several 
seconds  from  station  to  station  but  in  an  understandable  way.  Comparing 
the  records  at  ATL  from  events  3  and  9  one  can  see  that  the  Love  pulse 
is  faster  from  event  3.  From  the  map,  in  figure  2,  it  is  seen  that  the 
path  from  event  3  is  3/4  ocean  older  than  100  m.y.  whereas  the  path 
from  event  9  is  only  about  1/2  ocean  that  old.  Past  studies,  such  as 
that  by  Mitchel  &  Yu  (1980),  have  documented  the  growth  in  thickness  and 
velocity  of  the  lid  with  age  in  the  Pacific.  We  are  probably  seeing 
this  phenomenon  here.  Our  model  is  an  average  of  a  varying  lid  in  the 
older  Atlantic.  Event  3,  at  ATL,  has  a  faster  Love  wave  than  our  model 
predicts  indicating  the  lid  is  probably  even  thicker  than  100  km  in  the 
older  ocean,  while  the  paths  from  these  events  are  probably  through  a 
thinner,  lower  velocity  lid.  Never  the  less,  a  100  km,  4.75  km/sec  lid 
does  a  fair  job  over  these  paths,  as  it  did  for  the  data  at  nearin 
distances . 

In  figure  13  further  SS  data  is  presented  from  distances  beyond 
40°.  The  stations  used  for  this  profile  are  well  inland  but  again  the 
midpoints  are  within  old  ocean.  The  large  distances  used  here  and  the 
fact  that  the  stations  are  so  far  inland  eliminate  contamination  by  the 
Love  wave.  At  SHA  long  period  Love  wave  energy  is  apparent  but  it  is 
behind  the  SS  wave  and  does  not  present  a  problem.  Since  the  Love  wave 
is  not  important  here,  and  the  SS  phase  at  these  distances  is  dominated 
by  structure  well  below  the  lid,  we  computed  the  synthetics  in  this 
figure  using  the  more  efficient  WKBJ  technique.  No  Love  wave  is 


computed.  In  this  range  SS  passes  through  the  same  triplications  as  S 
from  20  to  26°.  The  branches  from  the  travel-time  curve  in  figure  7  are 
labeled  on  the  synthetics  and  the  data.  The  data  at  41  and  43°  are 
reproduced  many  times  from  other  events  and  at  other  stations.  At  41° 
two  arrivals  are  present,  the  first  is  from  the  400  km  discontinuity 
(branch  C),  the  second  from  the  660  km  discontinuity  (branch  E).  At  43° 
they  are  crossing  producing  a  large  simple  SS  wave.  Near  46°  the  back 
branch  of  the  400  km  discontinuity  begins  to  separate  out  and  can  be 
followed  to  53°.  The  660  back  branch  begins  to  separate  from  energy 
below  660  km  near  48°  and  it  too  can  be  followed  to  53°.  This  data  just 
confirms  the  previous  conclusions  but  also  illustrates  the  compatiblity 
of  S  and  SS  wave  modeling.  Many  regions  which  are  inacessible  to  an  S 
wave  upper  mantle  study  can  be  sampled  with  SS. 


Discussion 


The  model  derived  in  this  study  was  based  on  several  assumptions. 
The  lid  was  assumed  to  be  a  homogeneous  high  velocity  layer.  Structure 
in  the  lid  is  certainly  possible  but  our  technique  of  calculating 
synthetics  prohibited  using  more  than  a  couple  of  layers  due  to  the 
number  of  rays  needed.  We  have  concentrated  on  fitting  long  period 
data,  >  20  sec,  which  should  not  be  too  sensitive  to  the  fine  structure 


of  the  lid.  A  high  gradient  in  the  lid  would  necessitate  a  thicker  lid 
as  the  shadow  would  be  moved  to  nearer  distances.  Thus,  our  lid  is  just 
an  approximation  to  a  more  complicated  structure  which  could  be  thicker. 

The  velocities  just  below  the  lid  are  difficult  to  determine  due  to 
the  extensive  shadow  zone  caused  by  the  lid.  We  used  the  same  structure 
below  200  km  as  was  found  beneath  the  East  Pacific  Rise  in  paper  1,  but 
there  is  some  tradeoff  between  the  velocity  from  100  to  200  km  and  the 
velocity  from  200  to  400  km.  On  the  other  hand,  decreasing  the  velocity 
in  the  Atlantic  model  just  below  the  lid  to  values  found  beneath  the 
East  Pacific  Rise  would  allow  only  a  slight  increase  in  velocity  from 
200  to  400  km.  Models,  such  as  PREM  (Dzeiwonski  &  Anderson,  1981),  with 
large  velocity  jumps  near  200  km  and  corresponding  low  gradients  to 
400  km  can  not  satisfy  our  data.  The  back  branch  from  the  400  km 
discontinuity  is  too  fast  in  these  models.  Thus  we  feel  confident  in 
the  conclusion  that  the  shear  velocity  from  200  to  400  km  beneath  the 
old  Atlantic  is  far  closer  to  the  velocities  beneath  the  East  Pacific 
Rise  than  beneath  the  Canadian  shield. 

It  is  interesting  to  compare  our  results  to  other  studies  of  the 
upper  mantle.  Sipkin  and  Jordan  (1976,1980)  have  measured  ScS  and 
multiple  ScS  times,  attributing  their  variations  to  upper  mantle 
heterogeneity.  For  the  shield  model  (SNA),  with  a  35  km  thick  crust, 
ScS  is  4  sec  faster  than  that  predicted  by  the  Jef f reys-Bullen  model. 
The  rise  model,  TNA,  is  6.2  seconds  slower  and  the  Atlantic  model  is 
just  about  the  same  as  the  Jef f reys-Bullen  prediction.  For  the  two 
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oceanic  models  we  used  a  9  km  crust.  The  shield  model  fits  the  mean  of 
the  measurements  by  Sipkin  and  Jordan  (1976)  for  Precambrian  shields, 
though  they  have  a  2  sec  scatter  for  those  measurements.  Their  oceanic 
data  fall  between  -1  and  8  second  residuals.  These  are  mostly 
measurements  in  the  Pacific  ocean  with  varying  oceanic  age.  Our  two 
oceanic  models  bracket  their  data  fairly  well.  Thus  the  two  studies 
seem  compatible  as  far  as  generalizing  our  results  to  similar  tectonic 
provinces  on  other  parts  of  the  earth. 

Figure  14  shows  the  fundamental  Rayleigh  wave  phase  velocity 
dispersion  curves  for  the  three  models  derived  using  S  and  SS. 
Corrections  due  to  analestacity  were  included  assuming  that  our  bodywave 
models  are  appropriate  for  a  period  of  10  sec.  The  Q  structure  was 
taken  from  PREM  (Dzeiwonski  &  Anderson, 1981).  The  P  velocity  for  the 
shield  was  taken  from  Given  &  Helmberger  (1981)  and  for  the  East  Pacific 
Rise  from  Walck  (1983).  The  Atlantic  P  model  is  a  modified  version  of 
Walck's,  with  a  lid  corresponding  to  the  shear  wave  Atlantic  model  and 
the  same  velocities  below  200  km.  Again  a  crust  of  35  km  was  used  for 
the  shield  model  and  one  of  9  km  for  the  two  oceanic  models.  Data  from 
several  regional  surface  wave  studies  is  included  on  the*  curves.  The 
Weilandt  &  Knopoff  study  (1982)  is  for  a  path  along  the  East  Pacific 


Rise  and  agrees  very  well  with  the  dispersion  predicted  by  TNA  at  all 
periods  but  the  shortest.  The  longer  period  study  by  Kanamori  (1970)  is 
for  an  average  ocean  and  agrees  with  the  long  period  results  for  the  ATL 
model.  Finally,  the  Okal  study  (1978),  looked  at  paths  across  North 
America.  The  measured  phase  velocities  are  lower  than  those  predicted 


by  our  shield  model  but  the  Okal  paths  include  significant  portions  of 
the  slower  western  part  of  North  America.  The  long  period  studies  agree 
well  with  our  results.  There  are  4%  differences  in  phase  velocity  at 
150  sec  decreasing  to  1.5%  differences  at  270  sec  for  the  three  tectonic 
regions.  At  the  shorter  periods  (<120  sec)  the  agreement  with  other 
studies  is  worse.  In  particular,  the  results  of  Mitchel  &  Yu  (1980)  for 
old  ocean  are  much  lower  than  ours.  Similarly,  the  Brune  &  Dorman 
(1963)  study  of  the  Canadian  shield  found  lower  velocities  than  our 
shield  model  predicted.  This  could  be  due  to  anisotropy,  as  recently 
discussed  by  Anderson  &  Regan  (1983),  or  actual  lower  values  of  Q,  than 
we  used,  at  shallow  depths.  The  Mitchel  &  Yu  study  was  for  the  Pacific 
ocean  and  could  indicate  a  real  difference  between  the  old  Atlantic  and 
old  Pacific.  In  conclusion,  the  long  period  Rayleigh  wave  dispersion 
predicted  by  our  models  agrees  well  with  other  studies  in  other  areas  of 
the  earth,  indicating  good  agreement  with  our  deeper  structure.  The 
shorter  periods  show  large  variations  which  could  be  due  to  a  number  of 
causes. 


Conclusion 


We  have  modeled  the  waveforms  and  traveltimes  of  SH  wave  motion  in 
the  older  portions  of  the  northwest  Atlantic  to  investigate  the  upper 


mantle.  The  old  Atlantic  has  a  thick  high  velocity  lid  beneath  which 
the  structure  Is  very  similar  to  that  found  near  the  East  Pacific  Rise 
in  a  similar  study.  Using  the  same  technique,  the  Canadian  shield  was 
found  to  have  significantly  higher  velocity  than  both  oceanic  areas  to 
about  400  km  depth.  The  bottom  of  the  lid  does  not  appear  to  be  the 
depth  at  which  the  craton  and  oceans  become  indistinct. 

Below  400  km  the  model  derived  for  the  Canadian  shield  and  East 
Pacific  Rise  fit  the  Atlantic  data  quite  well.  Thus  we  feel 
heterogeneity  at  these  depths  is  very  small  throughout  North  America  and 
its  surrounding  oceans.  The  agreement  with  very  long  period  surface 
wave  studies  elsewhere  suggest  the  model  below  400  km  is  a  pretty  good 
average  for  the  earth. 
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Table  1 


Atlantic  Model 


Depth , km  Velocity, km/sec  Depth, km  Velocity, km/sec 


0.0 

3.7 

450.0 

5.09 

9.0 

3.7 

475.0 

5.135 

10.0 

4.75 

500.0 

5.19 

25.0 

4.75 

525.0 

5.24 

50.0 

4.75 

550.0 

5.29 

75.0 

4.75 

575.0 

5.345 

100.0 

4.75 

600.0 

5.395 

125.0 

4.54 

625.0 

5.445 

150.0 

4.49 

659.0 

5.50 

175.0 

4.47 

660.0 

5.91 

200.0 

4.46 

675.0 

5.98 

225.0 

4.46 

700.0 

6.05 

250.0 

4.48 

725.0 

6.13 

275.0 

4.51 

750.0 

6.20 

300.0 

4.57 

775.0 

6.22 

325.0 

4.63 

800.0 

6.24 

350.0 

4.68 

825.0 

6.26 

375.0 

4.73 

850.0 

6.275 

405.0 

4.78 

875.0 

6.29 

406.0 

5.00 

900.0 

6.305 

425.0 

5.05 

925.0 

6.32 

cMco^m^or^cooH 
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Figure  1.  Velocity  model  ATL  derived  In  this  study  compared  to 
profiles  of  the  Canadian  shield  (SNA)  and  the  tectonic  western  North 
America  (TNA). 

Figure  2.  Stations  and  events  used  In  the  study.  The  age  of  the 
oceanic  crust  In  the  region  Is  also  indicated. 

Figure  3.  The  step  responses  after  summing  various  sets  of  rays  In 
two  simple  structures.  The  final  panel  is  the  response  of  a  long-period 
WWSSN  Instrument  to  the  Green's  function  above  it. 

Figure  4.  Step  responses  at  12°  for  various  simple  models.  To  the 
left  are  the  structures  used,  in  the  middle,  the  step  response,  and  to 
the  right  the  long-period  WWSS  instrument  response  with  a  time  function 
and  attenuation.  The  dashed  line  indicates  the  same  reference  time. 

Figure  5.  Data  and  synthetics  for  two  oceanic  paths  at  nearin 
stations.  The  timing  is  absolute.  The  data  is  a  above  the  appropriate 
synthetic. 

Figure  6.  Data  and  synthetics  for  the  Bermuda  event  to  continental 
stations.  The  400  km  reflection  is  marked  by  an  arrow.  The  timing  is 
absolute. 

Figure  7.  Triplication  curve  for  the  Atlantic  model.  Dashed  lines 
represent  diffracted  or  tunnelled  energy. 

Figure  8.  A  short  profile  of  S  waves  in  the  old  Atlantic  with  a 
synthetic  profile  to  the  right.  The  triplication  branches  corresponding 
to  figure  7  are  indicated. 

Figure  9.  A  comparison  of  the  fits  of  the  Atlantic  model  to  data 
with  a  modified  Atlantic  model.  The  modified  model  has  shield 
velocities  from  200  km  to  400  km  and  this  changes  the  400  km  back-branch 
significantly  as  Indicated  by  the  arrows. 

Figure  10.  The  construction  of  synthetic  seismograms  from  26°  to 
40°.  The  first  column  shows  the  effect  of  just  the  S  rays,  the  middle 
has  crustal  and  lid  multiples  added  and  the  right  column  has  deeper  SS 
rays  added.  The  triplication  branches  are  in  reference  to  figure  7. 

Figure  11.  A  sensitivity  study  on  the  effect  of  lid  thickness  on 
the  front  of  the  fundamental  Love  wave.  Just  lid  and  crustal  multiples 
are  computed  for  various  lid  thicknesses. 

Figure  12.  Data  and  synthetic  comparison  from  26°  to  40°.  The 
synthetics  are  from  figure  10,  several  triplication  branches  from  figure 
7  are  indicated. 

Figure  13.  Data  and  synthetic  comparison  from  410  to  53Q  ^  gg 
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phase  is  passing  through  the  triplications  as  indicated  by  the 
travel-time  lines.  The  first  arrival  S  wave  is  a  useful  time  marker. 

Figure  14.  Fundamental  Rayleigh  wave  dispersion  curves  appropriate 
for  the  Atlantic  model  from  this  study,  the  shield,  and  tectonic  models 
found  in  paper  1  are  compared  with  data  obtained  from  the  literature. 
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Applications  of  the  Transmitted  Kirchhoff-Helmholtz  Method 
to  Transmitted  Body  Haves  and  Possible  Structural  Effects  at  NTS 


Patricia  Scott 
and 

Donald  V.  Helmberger 


ABSTRACT 


We  extend  the  Kirchhoff-Helmboltz  integral  method  to  calcu¬ 
late  acoustic  potentials  which  transmit  through  three  dimensional 
warped  boundaries.  We  specify  the  potentials  on  an  arbitrary  sur¬ 
face  with  Snell’s  law  and  plane-wave  transmission  coefficients  and 
numerically  integrate  their  contributions  at  a  receiver  via  the 
scalar  integral  representation  theorem.  The  method  is  appropri¬ 
ate  for  modeling  precritical  transmitted  potentials.  Results  from 
test  models  compare  well  with  optical  solutions  for  transmissions 
through  a  flat  interface.  We  model  the  effect  of  several  idealised 
crust-mantle  boundary  structures  on  teleseismic  P-waves  gen¬ 
erated  by  explosion  sources.  The  structures  are  all  upwarps  and 
are  designed  to  produce  travel  time  residuals  as  a  function  of  delta 
and  azimuth  which  have  the  same  magnitude  as  residuals  observed 
for  NTS  tests  within  Pahute  Mesa.  These  structures  consistently 
cause  complicated  low  amplitude  waveforms  which  arrive  early  and 
simple  high  amplitude  waveforms  which  arrive  late.  Thus  they 
cause  systematic  amplitude  variations  with  azimuth,  delta,  and 
source  location.  The  magnitude  of  this  variation  is  less  than  or 
equal  to  2}£.  This  factor  is  somewhat  less  than  the  observed  ab 
amplitude  variation  with  azimuth  of  Pahute  Mesa  tests;  however  it 
is  approximately  the  same  magnitude  as  the  observed  ab  variation 
at  a  given  station  as  a  function  of  test  location  within  the  mesa. 


I.  Introduction 


Despite  dramatic  improvements  in  the  level  of  sophistication  of  data 
analysis,  seismologists  still  cannot  deterministically  predict  many  observed 
amplitude  and  travel  time  anomalies  of  body  waves.  One  hypothesis  to  explain 
these  anomalies  is  the  presence  of  non-planar  velocity  discontinuities  such  as 
sedimentary  basins,  mountains,  and  faults  near  the  source  or  receiver.  In  order 
for  us  to  assess  the  importance  of  this  hypothesis,  we  need  a  technique  for 
predicting  the  impact  of  near  source  and/or  near  receiver  structural  complex¬ 
ity  on  far-fleld  waves.  In  this  paper,  we  present  such  a  method  based  on  the 
numerical  evaluation  of  the  Kirchhofl-Helmholtz  integral  with  use  of  modified 
tangent  plane  boundary  conditions.  This  method  calculates  the  response  of  a 
wave  which  is  transmitted  through  a  warped  boundary  between  two  acoustic 
media.  It  contrasts  from  Scott  and  Helmberger  (1983)  where  the  reflections 
from  a  warped  boundary  are  calculated.  In  this  paper  we  briefly  describe  the 
formalism  of  Kirchhofl-Helmholtz  method  for  the  transmitted  case.  Then,  as  an 
example  of  the  method,  we  model  the  observed  azimuthal  amplitude  and  travel 
time  anomalies  of  short  period  P  waves  from  NTS  blasts  as  a  result  of  a  geologic 
structure  at  the  Moho. 

II.  Formalism 

The  method  in  this  work  is  based  on  the  numerical  evaluation  of  the 
Kirchhofl-Helmholtz  integral  equation.  The  formalism  differs  slightly  from  that 
presented  in  Scott  and  Helmberger  (1983).  In  that  paper  we  discussed  the 
reflected  wave  solution;  however,  here,  we  state  the  Kirchhofl-Helmholtz  solution 
for  a  transmitted  potential.  We  also  qualitatively  discuss  the  assumptions 


involved  in  its  use. 


We  wish  to  calculate  a  transmitted  potential  at  pointy. in  a  homogeneous 


body  Vg  resulting  from  an  incident  source  potential  located  at  point  x,  in  a ' 


homogeneous  body  K,.  The  boundary  between  the  two  bodies  is  dV.  The  sound 


speeds  and  densities  of  Vx  and  Vg  are  ax  and  ag  and  px  and  pg  respectively. 


From  the  scalar  integral  representation  theorem,  we  write  the  solution  for  <pg  at 


a  pointy. off  the  boundary  dV,  within  Vg,  and  at  a  time  t  as 


w«CfcO=/[c*(i )-*(*.« )*IftCExO]*ii8(aL)«<5  .  d) 

dv 


Here  •  denotes  convolution  and  •  denotes  a  vector  dot  product.  Gg  is  the  funda¬ 


mental  singular  solution  of  the  scalar  wave  equation. 


if  y2£?e = 6{t 


In  addition  nz  is  the  outward  pointing  normal  of  Vg.  $  and  are  the  poten¬ 


tial  and  its  normal  derivative  on  the  surface  dV  in  Vg.  We  display  the  geometry 


in  Figure  1(a)  for  this  problem. 


Equation  (1)  is  exact  for  the  initial  conditions  that  <pg(j£j 0)  and  ^(j^O)  equal 


zero  throughout  V2.  The  derivation  can  be  found  in  Mao  and  Pao  (1971)  or  Strat¬ 


ton  (1941).  To  obtain  equation  (1),  one  requires  that  dV  is  a  reasonably  smooth 


surface. 


We  now  specify  Gg  and  $.  For  a  homogeneous  medium,  it  is  sufficient  to  use 


G  /,  ,  t  \  -  .  _  _  jg-Lj 
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If  the  incident  field  in  V,  results  from  an  isotropic  point  source  at  £ then  $  is 


approximated  by 


,  Tf(t-  rt) 

T>  = 


\z-Zo 


Here  T  is  the  acoustic  plane  wave  transmission  coefficient  for  a  flat  interface  and 
is  dependent  on  the  local  incidence  angle  at  each  point,  /(f)  is  the  time  func¬ 
tion  of  the  incident  source  function.  The  function  in  equation  (4)  approximates 
♦  well  if  the  incident  source  field  is  of  sufficiently  high  frequency  such  that  every 
point  on  the  surface  transmits  the  incident  pulse  as  though  there  were  an 
infinite  plane  tangent  to  the  surface  at  that  point.  Then  the  amplitude  and  the 
phase  on  the  surface  can  be  described  locally  by  plane  wave  transmission 
coefficients  and  Snell's  law.  The  value  of  the  potential  at  one  point  is  indepen¬ 
dent  of  the  values  at  other  points.  Hence  the  contributions  of  diffractions  and 
multiple  scattering  to  the  potential  are  neglected. 

d  Cr% 

We  now  estimate  the  normal  derivatives  — -  and  — - . 


-TJ(t-rx)  dr. 


a£j_  s  -6(t  -Tg)  8t2 

fln2  4  rr  Le.~JL  I  dn2 

The  dot  over  the  functions  in  (5)  and  (B)  signifies  time  derivative.  We  approxi¬ 
mate  the  normal  derivatives  by  assuming  the  amplitudes  of  Gz  and  $  vary  slowly 
on  the  surface  relative  to  the  phase.  Hence  we  can  discard  the  terms 

T~6(t  -Tg)-~ — (1 — - — r-)  and  f(t- t,)-^ — (1 — - — 7-).  Substitution  of  equations 
4rr  ong  lx.— il  an 2  \j£.-£o  I 


(3).  (4),  (5)  and  (6)  into  equation  (l)  yields 


=  1  r  r/(f~T»-T2)  (6t2  6ti 

'  ®  4tt  »  r  r  '  An  -  - 


where  r  —  |x— xl.  the  distance  from  the  surface  to  the  receiver,  and  r0  =  |ar— |, 
the  distance  from  the  source  to  the  surface.  We  note  that  the  discarded  parts  of 

normal  derivatives  are  proportional  to  — — g-and  — — g-.  For  problems  computed 

in  this  study,  the  distance  from  the  source  to  the  surface  averages  40  kilome¬ 
ters  and  the  distance  from  the  surface  to  the  receiver  averages  20,000  kilome¬ 
ters.  The  — — j-  and  ^  terms  are  6.25xl0“18  and  3.1x10“®,  respectively,  and 

are  small  relative  to  the  term  — — in  equation  (7),  which  is  1.25xl0“7.  We  there- 

t0t 

fore  discard  these  terms  with  confidence. 

<Jt,  dr 9 

To  estimate  — - and  - - ,  we  recall  that  the  gradient  of  the  phase  is 

an2  dnz 

parallel  to  the  normal  of  the  wavefront  and  has  units  of  slowness.  Thus, 


cos0, 


where  cos0i  is  the  cosine  of  the  angle  between  the  normal  to  the  refracted  wave- 
front  and  the  normal  to  the  surface  n2.  It  is  equal  to 


cosB j  =  (1  -  — 5—  sin2  0i )& 

of 


where  0t  is  the  local  incident  angle  calculated  by 


cosO*  = 


(x-Xt) 

•  • 

r , 


Hz  ■ 


(ID 


cos02  is  the  cosine  of  the  angle  between  the  normal  n2  and  a  ray  connecting  the 
surface  and  the  receiver.  Thus 


cos02 


(12) 


Figure  1(b)  shows  a  detailed  picture  of  these  angles.  Substitution  of  the  cosine 
factors  (10)  and  (12)  yields 


Vz 


COS02 


1  r  •  ,  .. COS 


COS0] 


r.rtt! 


■)  dS  . 


(13) 


The  method  for  the  calculation  of  equation  (13)  is  discussed  in  Scott  and 
Helmberger  (1983).  We  simply  calculate  the  integral  as  a  summation  of  single 
point  evaluations  of  the  integrand.  This  method  of  integration  requires  that  the 
elements  which  comprise  the  surface  be  small  in  length  compared  to  the 
incident  source  wavelength.  As  in  Scott  and  Helmberger  (1983)  we  obtain  the 
numerical  ramp  response  from  this  integration  and  convolve  it  with  the  analyti¬ 
cal  third  derivative  of  a  Haskell  isotropic  source.  Thus  we  obtain  £2,  the  time 
derivative  of  the  potential. 

These  calculations  are  appropriate  for  precritical  transmissions  in  a  linear 
acoustic  media.  We  do  not  allow  the  transmission  coefficient  to  be  complex. 
When  the  incident  angle  0*  exceeds  the  critical  angle,  $  equals  zero.  This  boun¬ 
dary  value  is  consistent  with  geometric  ray  theory,  but  is  not  a  realistic 
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shadowing  function.  Rather  we  expect  $  and 


8$ 

— - to  vary  smoothly  across  the 

OTLg 


shadow  boundary  if  the  surface  is  reasonably  smooth.  However,  we  argue  that, 
under  the  circumstances  considered  here,  the  postcritical  incidence  portions  of 
the  integral  contribute  to  the  summed  response  of  the  transmitted  potential  at 
times  much  later  than  the  specular  portions  and  .therefore,  are  unimportant. 


As  a  check  of  the  method,  we  compute  the  transmitted  response  of  a  wave 
propagating  through  a  planar  boundary  separating  two  volumes  of  different 
sound  speeds.  We  then  compare  the  numerical  maximum  amplitude  with  that 
obtained  from  the  following  first  motion  formula: 


<Pt  =  TLf(t 

al  <*2 


(14) 


Here  L  is  the  spreading  coefficient  (Langston,  1977;  Hong.  1978) 


L  = 


!(rJ-+  zr~)( 


Vv 


VvZ  vfTJv, 


22*7*1  v 

'  o  «  / 

v2  *7v2 


(15) 


where  and  z2  are  the  vertical  distances  of  the  source  and  receiver  respec 
tively  from  the  boundary,  T)^  is 


*7  vi  = 


(16) 


where  p,  is  the  ray  parameter. 

Figure  2  shows  examples  of  this  comparison.  We  have  computed  transmit¬ 
ted  potentials  for  an  incident  isotropic  source  which  is  the  first  derivative  of  a 
Haskell  source  with  parameters  (B=2,  K=10).  The  velocity  and  density  model 
used  for  the  comparison  calcule '.ions  are  shown  in  the  top  of  Figure  2(b). 
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We  show  two  KirchhofI  synthetics  in  Figure  2(a)  to  demonstrate  the  nature 
of  truncation  phases  which  can  contaminate  the  synthetics.  These  phases  arrive 
approximately  3$  seconds  after  the  first  arrival  in  both  synthetic  A  and  syn¬ 
thetic  B  and  they  are  artifacts  of  the  technique.  In  synthetic  A.  the  phase  is  a 
result  of  the  finiteness  of  the  grid.  A  diagram  to  the  left  of  the  synthetic  shows 
this  effect.  The  grid  is  a  square  with  a  length  of  150  kilometers.  The  source  is 
500  kilometers  above  the  cent  .r  of  t;  a  grid  and  the  receiver  is  1000  kilometers 
below  the  source.  From  this  diagram,  we  observe  that  the  edge  interferes  at  a 

time  t2-  •  The  geometric  ray  arrives  at  a  time  /,=  ^-+  — —  .  Hence 
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the  truncation  phase  arrives  3.6  seconds  later  than  the  first  arrival  in  synthetic 

A. 


In  synthetic  B,  the  phase  is  a  result  of  the  shadowing  function.  We  use  the 
same  grid  to  calculate  synthetic  B  as  for  synthetic  A;  however  the  source  is  167 
kilometers  above  the  interface.  For  this  velocity  model,  the  local  angle  between 
the  incident  ray  and  the  normal  to  the  surface  exceeds  critical  when  the  dis¬ 
tance  from  the  center  of  the  grid  exceeds  50  kilometers.  $  on  the  surface  is  0 
beyond  this  distance.  This  abrupt  change  in  boundary  conditions  introduces  a 
truncation  phase  into  the  synthetic.  From  the  diagram  to  the  left  of  synthetic, 


'  t- 

,',v 

"l- 

£ 
& 

I 


& 

$ 


we  see  that  this  phase  arrives  4  seconds  later  than  the  geometric  arrival.  The 
truncation  phase  in  synthetic  A  caused  by  grid  finiteness  does  not  constitute  a 
problem.  If  it  contaminates  the  phase  of  interest,  we  can  enlarge  the  grid 
appropriately.  However,  the  truncation  phase  in  synthetic  B  caused  by  the 
boundary  conditions  fundamentally  restricts  the  source-receiver  geometries  we 
can  investigate. 
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Figure  2(b)  shows  a  profile  of  KirchhofT  synthetics  for  a  source  500  kilome¬ 
ters  above  the  interface  and  five  receivers  500  kilometers  below  the  interface. 
The  horizontal  distance,  x,  of  the  receivers  ranges  from  0  kilometers,  directly 
underneath  the  source,  to  755  kilometers.  The  two  columns  next  to  the  synthet¬ 
ics  contain  the  numerical  peak  amplitudes  and  the  predicted  amplitude  from 
equation  (15).  The  agreement  is  good.  We  cannot  calculate  a  response  past  x  = 
755  kilometers  because  a  truncation  phase  resulting  from  the  boundary  condi¬ 
tions  on  the  interface  starts  to  interfere  with  the  direct  arrival.  We  must  always 
take  care  to  avoid  such  contamination. 

III.  NTS  Structure  (An  Example  of  Near-Source  Effects) 

We  now  apply  the  method  by  modeling  the  effects  of  idealised  Moho  struc¬ 
tures  on  transmitted  teleseismic  P  waves  generated  by  nuclear  tests  in  Pahute 
Mesa,  Nevada  Test  Site.  We  wish  to  ascertain  whether  focusing-defocusing  by 
structure  on  the  Moho  explains  the  unusual  behavior  of  amplitudes  from  these 
tests. 

We  review  these  anomalous  observations  of  short  period  P  waves  from 
Pahute  Mesa.  Figure  3  is  a  plot  of  1200  ab  amplitude  measv  ements  from  25 
tests  within  Pahute  Mesa  as  a  function  of  station  location  from  Lay  et  al.. 
(1983a).  The  ab  amplitudes  are  measured  from  the  first  peak  to  the  first  trough. 
They  are  corrected  for  geometric  spreading,  the  instrument  gain  at  1  second 
and  event  size,  following  a  procedure  developed  by  Butler  (1984).  The  ampli¬ 
tudes  are  relative  to  a  master  event  which  minimises  the  overall  scatter  of  the 
data. 
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The  data  has  two  important  features.  First,  the  relative  amplitudes  range 
from  .13  at  station  TRI  to  5.1  at  station  SHK.  This  variation  is  nearly  a  factor  of 
40.  Most  stations  between  the  azimuths  0°  and  60®  have  significantly  lower 
amplitudes  than  those  between  60®  and  120®.  Secondly,  the  relative  amplitudes 
at  a  given  station  vary  by  a  factor  of  2}£  as  a  function  of  event  location  within  the 
mesa.  The  latter  variation  clearly  originates  from  a  near*  source  mechanism 
because  the  events  are  separated  by.  at  most,  15  kilometers. 

If  one  calculates  the  mean  relative  amplitude  at  each  station,  then  the 
overall  amplitude  variation  with  azimuth  reduces  to  a  factor  of  12  (Lay  et  al., 
1983a).  The  next  two  figures  suggest  that  this  large  amplitude  scatter  is  also 
caused  by  a  near-source  mechanism.  Figure  4,  from  Lay  et  al.,  (1983a),  shows 
the  azimuthal  pattern  of  relative  amplitudes  for  GREELEY,  an  event  within  the 
mesa,  and  FAULTLESS,  an  event  100  kilometers  outside  the  mesa.  Although  both 
events  have  comparable  yields,  their  azimuthal  patterns  differ  substantially. 
This  difference  is  particularly  obvious  between  0®  and  90®.  Figure  5  displays 
plots  from  Lay  et  al.,  (1983b)  which  enhance  the  difference  between  patterns  of 
events  in  the  mesa  and  events  outside  the  mesa.  These  plots  are  ratios  of  ampli¬ 
tudes  of  three  events  outside  the  mesa  (FAULTLESS,  PILEDRIVER.  and  B1LBY) 
divided  by  the  average  mesa  amplitudes.  These  ratios  are  an  approximate 
measure  of  a  near-source  anomaly  if  the  FAULTLESS,  PILEDRIVER,  and  BILBY 
patterns  are  only  influenced  by  path  and  receiver  effects  and  are,  therefore, 
constant  as  a  function  of  azimuth.  Furthermore,  the  path  and  receiver  effects 
must  be  characterised  by  multiplicative  factors.  Because  the  ratio  patterns  for 
all  three  events  are  similar,  these  assumptions  are  probably  true.  Therefore  the 
factor  of  13  variation  of  these  ratios  between  0®  and  120®  is  roughly  an  estimate 
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of  the  near-source  anomaly  at  the  mesa. 

To  see  if  this  amplitude  variation  correlates  with  waveform  changes,  we  plot 
in  Figure  6  several  seismograms  at  stations  between  30s  and  100°  which 
recorded  both  FAULTLESS  and  GREELEY.  The  top  and  bottom  seismograms  are 
recordings  of  FAULTLESS  and  GREELEY,  respectively,  with  their  absolute  ab 
amplitudes  in  millicrons,  corrected  for  instrument  gain  only.  There  is  no  obvi¬ 
ous  waveform  differences  in  the  GREELEY  records  which  correlate  with  the 
dramatic  ab  amplitude  changes.  Furthermore  we  do  not  see  any  obvious 
difference  in  frequency  content  and/or  complexity  between  low  stations  and 
high  stations  for  either  event.  However  there  are  some  systematic  differences 
between  GREELEY  and  FAULTLESS  seismograms.  A  shoulder  occurs  2  to  3 
seconds  after  the  first  arrival  on  GREELEY  records  (e.g.  STU,  PTO,  MAL,  STJ,  OTT, 
GEO,  and  ATL).  Lay  has  also  seen  these  arrivals  for  other  mesa  events  (Lay  et 
al. ,1983b).  No  such  arrival  is  apparent  on  the  FAULTLESS  seismograms.  Also  the 
width  of  the  first  pulse  of  GREELEY  seismograms  is  narrower  than  those  of 
FAULTLESS  seismograms  at  a  few  stations  (e.g.  SJG,  ATL,  BLA.  GEO,  SCP,  and 
STU).  Both  phenomena,  though,  occur  throughout  the  azimuthal  range  and  do 
not  correlate  with  the  ab  amplitude  changes. 

The  data  demonstrates  that  near-source  anomalies  cause  a  variation  of  2% 
of  relative  ab  amplitudes  at  a  given  station  as  a  function  of  event  location  within 
the  mesa.  Moreover,  near-source  anomalies  also  cause  part  of  the  large  ab 
amplitude  variation  with  azimuth  (or  station  location)  from  mesa  events.  We 
cannot  completely  eliminate  contamination  of  the  azimuthal  pattern  by  path 
and  near-receiver  effects.  Certainly  near-receiver  effects  can  be  as  large  as 
those  observed  for  the  Pahute  mesa  tests  (Butler,  1984).  Yet  the  similarity  of 
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the  ratio  patterns  of  FAULTLESS,  BILBY,  and  P1LEDRIVER  suggests  that  the  pat¬ 
tern  for  mesa  tests,  seen  in  Figure  3,  is  dominated  by  a  near-source  mechanism. 
Finally  the  variation  of  relative  ab  amplitudes  with  azimuth  does  not  correlate 
with  any  obvious  waveform  changes  for  a  typical  mesa  event,  GREELEY.  There  is 
no  definitive  evidence  to  determine  whether  ab  variations  correlate  with  travel¬ 
time  residuals. 

In  this  paper  we  assume  that  all  the  observed  amplitude  anomalies  result 
from  near-source  mechanisms.  We  then  test  the  hypothesis  that  structure  on 
the  Moho,  consistent  with  travel-time  residuals,  focuses  or  defocuses  P-waves 
enough  to  produce  the  magnitude  of  the  amplitude  anomaly.  There  are  alterna¬ 
tive  near-source  explanations  for  these  anomalies.  In  addition,  to  the  focusing- 
defocusing  hypothesis,  workers  {Lay,  et  al.,  1983;  Wallace  et  al..  1983  )  postulate 
that  the  movement  of  faults  associated  with  nuclear  blasts  causes  a  superposi¬ 
tion  of  distributed  or  point  double-couple  sources  with  the  isotropic  bomb 
source.  The  amplitude  anomalies  are,  then,  the  radiation  pattern  caused  by  a 
double-couple  source.  Longer  period  studies  of  Love/Rayleigh  ratios.  Pnl,  P  and 
S  waves  (Aki  emd  Tsai,  1972;  Wallace  et  al.,  1983;  Nuttli,  1969)  generated  by  these 
blasts  support  the  latter  hypothesis.  However,  we  speculate  that,  as  the  fre¬ 
quency  content  of  the  signal  increases,  the  role  of  lateral  near-source  structure 
in  distorting  amplitudes  becomes  more  important.  From  travel  time  residual 
studies  (Minster  et  al.,  1981;  Spence,  1974)  workers  have  deduced  that  there  is  a 
high  velocity  zone  directly  beneath  the  Silent  Canyon  Caldera  in  the  mesa  which 
extends  down  to  100  kilometers.  Such  a  velocity  structure  may  cause  ampli¬ 
tudes  which  deviate  from  those  predicted  by  a  spherically  symmetric  Earth 
model. 
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To  investigate  how  geology  can  aflect  amplitudes,  we  presume  that  the 
apparent  velocity  variations  deduced  from  the  travel  time  residuals  are  a  man¬ 
ifestation  of  Moho  topography.  We  exclude  from  consideration  the  impact  of  the 
Silent  Canyon  volcanics  on  transmitted  P  waves  because  both  Spence  and  Min¬ 
ster  correct  the  residuals  statically  for  these  low  velocity  rocks;  thus  the  resi¬ 
dual  patterns  are  not  a  result  of  the  caldera.  In  any  case,  we  cannot  readily 
model  a  feature  so  close  to  the  source.  If  we  place  a  strong  velocity  discon¬ 
tinuity  .such  as  that  between  volcanic  and  granitic  rocks,  closer  than  10  kilome¬ 
ters  to  the  source,  we  generate  a  truncation  phase  which  interferes  with  the 
transmitted  P  phase. 

We  describe  the  Earth  with  a  two  layer  crust-mantle  velocity  model.  Tin 
velocity  of  the  upper  layer  is  6.5  km/sec  and  that  of  the  lower  layer  is  B  km/sec. 
The  depth  of  the  interface  is  45  kilometers.  The  receivers  are  located  at  dis¬ 
tances  such  that  the  —•  amplitude  decay  corresponds  to  spreading  at  telese- 

R 

ismic  distances  between  60°  and  70®  for  a  JB  Earth  (Langston  and  Kelmberger, 
1975). 

The  number  of  ways  to  distort  the  Moho  is  infinite.  We,  therefore,  restrict 
ourselves  to  a  few  three  dimensional  topographies  where  the  maximum  height  of 
the  anomaly  is  10  kilometers  and  the  maximum  width  is  approximately  25 
kilometers.  The  choice  of  these  values  is  based  on  both  the  Minster  et  al.  (1991) 
and  Spence  (1974)  studies.  They  find  an  advance  of  *.25-  4  seconds  for  nearly 
vertical  rays  from  shots  within  the  caldera.  As  these  rays  become  shallow,  this 
advance  lessens  or  disappears  completely.  From  crude  calculations,  we  esti¬ 
mate  that  10  kilometers  of  upward  relief  on  the  Moho  will  produce  the  required 
timing  anomalies  of  these  fays.  Furthermore,  we  confine  the  relief  laterally  so 


that  rays  exiting  the  mesa  at  shallow  angles  are  unaffected  by  the  structure.  We 
recognize  that  these  structures  are  extreme.  However,  if  we  cannot  produce  the 
observed  amplitude  anomalies  with  these  topographies,  we  can  discard  struc¬ 
ture  on  the  Moho  as  the  dominant  cause  of  these  anomalies. 

Of  the  infinite  number  of  structures,  we  arbitrarily  select  four  examples 
with  height  c  =  10  kilometers  and  width  w  =  25  kilometers.  These  topographies 
are  described  by  simple  analytical  formulas  and  are  convenient  to  use.  The 
topographies  with  their  labels  are  as  follows: 

Z  =  Zc9n  +  |~(cos(2rr((r-ii;/2)/ti;))-l)  if 

Upwarp  :  (17) 

Z  -  Zeon  if  r>f- 

— 4  0C3 r/ 

Exponential  :Z  =  Zem  ~  ce  2  (17a) 

Z  =  Z„n  -c  ifr«^- 

Plug  :  Z  =  Zcon  +  |<cos(2rr((r-|--  S)/w))-l)  if  f-<r s  |-  +  5  (17b) 

Z  ~~  Zem  if  r>-^— +  5 

Sine  :  Z  =  Z„n  ~  c  sine  [  21  ~^91r  ]  .  (17c) 

Here  r  =  V(x  -xc  )z  +  (y-yc  )z  *s  the  horizontal  distance  of  each  point  on  the 
surface  from  the  midpoints  of  the  grid  ( zc,ye ).  Zeon  is  the  baseline  level  of  the 
Moho  and  is  45  kilometers  for  all  the  calculations.  The  values  of  the  constants  in 
the  exponential  and  the  sine  bumps  confine  the  anomaly’s  width  to  approxi¬ 
mately  25  kilometers.  A  schematic  cross  section  of  each  topography  is  shown  in 
Figure  7.  All  the  structures  are  symmetric  with  azimuth. 
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Initially  an  isotropic  source  is  directly  above  the  center  of  the  structure; 
thus,  the  transmitted  potential  is  only  a  function  of  x  and  z.  The  source  is  45 
kilometers  above  the  baseline  of  the  Moho.  The  transmitted  potential  is  calcu¬ 
lated  at  receivers  which  are  20,000  kilometers  below  the  source.  The  horizontal 
distance  of  the  receivers  from  the  center  of  the  topography  ranges  from  0  to 
7,000  kilometers.  Figure  7  shows  the  transmitted  potential  and  the  peak  ampli¬ 
tude  as  a  function  of  x  in  increments  of  1,000  kilometers  for  each  of  the  four 
topographies.  In  addition,  the  responses  for  a  wave  which  transmits  through  a 
planar  boundary  are  displayed  in  the  first  column.  The  corresponding  take-off 
angle  for  the  flat  boundary  synthetics  are  to  the  left  of  the  column.  By  compar¬ 
ing  these  synthetics  with  those  in  the  other  columns,  we  can  determine  how 
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much  distortion  of  the  waveform  is  caused  by  each  structure.  The  synthetics  in 
Figures  7  and  10  do  not  include  a  Q  and  instrument  operator  or  a  reflected  pP 
phase.  Although  these  effects  are  important,  we  want  to  examine  amplitude  and 
waveform  distortions  caused  by  structure  with  a  simple  input  pulse.  The  ringing 
caused  by  an  instrument  or  pP  may  mask  the  presence  of  multiple  arrivals 
caused  by  the  topography. 

All  the  structures  cause  intriguing  changes  in  the  waveforms  and  arrival 
times  of  the  synthetics.  The  waveform  features  originate  from  timing  changes 
caused  by  each  topography.  Each  point  of  the  Kirchhoff  synthetic  originates 
from  elements  which  are  illuminated  by  the  source  and,  in  turn,  illuminate  the 
receiver  at  a  total  travel  time,  <j=t,+t8.  We  can  associate,  with  each  element  of 
the  grid,  a  value  of  total  travel  time,  o(x,y).  The  total  travel  time  function  on 
the  surface  depends  on  the  source  location,  the  receiver  location,  and  the  sur¬ 
face  geometry.  Figure  B  shows  examples  of  this  function.  Here  we  calculate 
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o(x,y)  for  elements  •which  make  up  a  flat  interface  8(a),  an  upwarped  interface 
B(b).  a  plug  8(c),  and  a  sine  function  B(d).  In  each  example  the  source  is  directly 
above  the  center  of  the  structure  and  is  45  kilometers  above  the  baseline  of  the 
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interface.  The  receiver  is  20,000  kilometers  directly  below  the  source.  The  con¬ 
tours  of  constant  total  travel  time  are  projected  onto  the  topography  (top 
figure  in  B(a).(b),(c),  and  (d)).  We  also  display  these  contours  as  a  function  of  x 
and  y  (bottom  figure  in  8(a),(b),(c),  and  (d)).  The  contours  are  circles  because 
of  the  particular  source-receiver  geometry.  For  the  sake  of  brevity,  we  only 
show  that  portion  of  the  grid  which  contributes  to  the  initial  second  of  each  Kir- 
chhoff  synthetic.  The  synthetics  which  correspond  to  these  total  travel  time 
functions  are  also  shown  (middle  figure  in  8(a), (b),(c), and  (d)).  These  figures 
show  how  structure  on  the  interface  distorts  the  total  travel  time  contours  and, 
as  a  result,  produces  multiple  arrivals  in  the  synthetic. 

We  examine  this  effect  in  detail.  The  contours  are  in  intervals  of  .125 
seconds  as  are  the  tick  marks  below  the  synthetics  in  Figure  8.  The  geometric 
arrival  time  occurs  at  the  center  of  the  contour  plot.  Thus,  by  counting  the  con¬ 
tours,  we  can  estimate  the  cumulative  area  of  the  surface  which  contributes  to 
the  synthetic  at  a  given  time.  We  deduce,  from  Figure  8.  that,  approximately, 

fPL  <ie> 

A(t)  is  the  amplitude  of  the  response  at  time  t.  S  is  the  total  area  of  the  surface 
which  contributes  to  the  response  at  time  t.  For  example,  the  initial  .375 
seconds  of  the  synthetic  from  a  flat  interface  results  from  a  rapid  increase  in 
the  cumulative  area  of  the  surface  which  is  illuminated  between  t=.125  and 
t=.25  seconds.  After  t=.25  seconds,  the  area  of  the  surface  is  illuminated  at  a 


constant  rate.  Thus,  the  resultant  synthetic  can  be  viewed  as  a  convolution  of 
the  source  time  function  with  a  step  function  which  starts  between  t=.125  and 
.25  seconds. 


We  quantify  this  statement  by  following  an  approach  developed  by  Hilter- 
man  (1975)  and  Haddon  and  Buchen  (1961).  The  symmetry  of  the  source- 
receiver  geometry  and  the  surface  geometry  allows  us  to  recast  the  integral 
(13)  as  a  one  dimensional  integration  with  respect  to  total  travel  time,  a.  If  the 
transmission  coefficient  varies  slowly  over  the  surface,  then 
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flj  and  f)2  are  the  time  derivatives  of  modified  solid  angles.  Q,  is  a  modified  solid 
angle  with  vertex  at  the  source  subtended  by  the  surface  S,  and  fl2  as  a  modified 
solid  angle  with  vertex  at  the  receiver  subtended  by  S. 

We  now  examine  the  origins  of  the  multiple  arrivals  in  the  potentials  which 
propagate  through  interfaces  with  structure.  For  example,  the  amplitude  and 
frequency  content  of  synthetics  from  the  upwarp  (column  3  of  Figure  7)  are 
controlled  by  the  interference  of  two  pulses.  The  travel  time  contours  in  Figure 
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8(b)  for  t.he  upwarp  differ  considerably  from  those  of  a  flat  interface  8(a).  Far 
less  of  the  upwarped  surface  is  illuminated  within  .25  seconds  of  the  geometric 
arrival  time.  Furthermore,  the  upwarp  topography  causes  subtle  changes  of  the 
width  between  travel  time  contours.  There  axe  two  locations  where  this  change 
occurs:  1)  at  the  top  of  the  upwarp  and  2)  at  the  edge  of  the  upwarp.  The  first 
pulse  in  this  synthetic  originates  from  the  elements  in  the  first  location  while 
the  second  pulse  originates  from  the  second  location.  Because  the  ring  of  ele¬ 
ments  w-hich  contribute  to  the  second  pulse  has  a  larger  area  than  that  of  ele¬ 
ments  which  contribute  to  the  first  pulse,  the  second  pulse  is  larger  than  the 
first  pulse. 

As  the  receivers  move  away  from  the  center,  the  maximum  amplitudes 
decrease  as  a  result  of  the  interference  of  the  two  pulses.  We  destroy  the  sym¬ 
metry  of  the  surface  illumination  by  moving  the  receivers  horizontally.  The 
illumination  of  elements,  which  initially  was  simultaneous,  now  occurs  at  slightly 
different  times  and  causes  destructive  interference;  this  destructive  interfer¬ 
ence  causes  a  reduction  in  peak  amplitudes  and  the  broadening  of  the  pulse 
widths  for  both  phases.  Moreover,  as  the  receiver  moves  out  laterally,  the 
planar  part  of  the  boundary  becomes  more  important  in  controlling  the  ampli¬ 
tude  of  the  transmitted  pulse.  Hence  the  amplitudes,  travel  times,  and  the 
waveforms  of  distorted  pulses  approaches  those  of  a  pulse  which  has  propagated 
through  a  flat  boundary.  This  phenomena  is  present  in  almost  all  the  top  syn¬ 
thetics  in  Figure  7. 

Other  intriguing  features  are  present  in  the  synthetics  shown  in  Figure  7. 
The  potentials  transmitted  through  the  sine  and  exponential  bumps  shown  in 
columns  3  and  5  have  an  apparent  delay  which  is  not  seen  in  the  other 


synthetics.  These  two  topographies  drop  in  height  near  the  peak  more  rapidly 
than  does  the  upwarp  topography.  Consequently,  fewer  elements  are 
illuminated  and  contribute  to  the  transmitted  potential  at  times  near  the 
geometric  arrival  time.  This  is  illustrated  for  the  sine  topography  in  Figure  B(d). 
Hence,  the  amplitude  near  the  geometric  arrival  time  is  lower  than  amplitudes 
at  later  times. 

The  opposite  is  true  for  the  synthetics  of  waves  which  are  transmitted 
through  a  plug.  They  are  displayed  in  column  4  of  Figure  7.  The  topography  and 
travel  time  contours  for  the  bottom  synthetic  from  this  column  are  shown  in 
Figure  B(c).  This  figure  shows  that  more  elements  are  illuminated  and  contri¬ 
bute  to  the  response  near  the  geometric  arrival  time  for  this  topography  than 
for  the  upwarp.  sine,  and  exponential  topographies.  The  resultant  synthetic  is 
made  up  of  two  pulses  of  equal  size.  Each  pulse  has  the  amplitude  and  shape  of  a 
wave  which  has  transmitted  through  a  planar  interface.  The  plug  is  essentially 
comprised  of  two  planar  interfaces,  one  at  z  =  45  kilometers  and  the  other  at  z 
=  35  kilometers.  The  edges  of  the  plug  have  been  tapered  to  avoid  a  shadowing 
problem.  While  the  difference  in  the  interface  depths  does  not  alter  the  ampli¬ 
tudes  of  the  pulses,  it  does  change  the  arrival  times.  This  slight  separation  in 
arrival  time  causes  the  observed  interference  pattern.  As  the  receiver  moves 


out  laterally,  the  timing  between  the  two  pulses  changes.  The  pulse  width  of  the 
first  arrival  narrows  while  the  width  of  the  second  one  broadens.  By  transmit¬ 
ting  a  wave  through  such  a  structure,  we  vary  the  maximum  amplitude  of  the 
synthetics  by  a  factor  of  2%. 

We  confirm  that  these  structures  are  approximately  producing  the  correct 
travel  time  anomalies.  We  plot  the  residuals,  in  addition  to  the  peak  amplitudes 
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of  the  synthetics,  as  a  function  of  distance  to  discern  any  systematic  relation¬ 
ship  between  the  two  parameters.  We  also  plot  the  amplitude  of  the  first  pulse  if 
the  synthetic  is  made  up  of  multiple  arrivals.  This  amplitude  is  measured  from 
the  start  of  the  synthetic  to  the  first  peak.  The  plots  are  displayed  in  Figure  9. 
The  travel  time  "residuals"  are  defined  by  the  difference  between  the  arrival 
times  of  transmissions  through  a  bumpy  surface  and  the  times  of  transmissions 
through  a  flat  surface.  Where  there  is  an  apparent  delay  in  the  synthetics  such 
as  in  those  from  the  exponential]  and  sine  bumps,  we  measure  the  arrival  time  at 
the  start  of  the  upswing.  The  amplitudes  are  uncorrected  for  geometric  spread¬ 
ing.  The  change  of  amplitude  from  spreading,  seen  in  the  synthetics  in  column  1, 
is  negligible  in  the  distance  range  of  interest;  hence  no  correction  is  necessary 

The  "residuals"  in  Figure  9  produced  by  these  structures  behave  in  a 
predictable  fashion.  The  transmitted  potentials  which  propagate  vertically  to 
stations  between  0  to  4,000  kilometers  experience  the  most  advance.  The 
exceptions  to  this  behavior  are  residuals  from  the  sine  and  exponential  synthet¬ 
ics.  We  know  there  is  some  energy  arriving  at  these  nearly  vertical  stations  with 
.3  second  advance  from  the  previous  discussion.  However,  because  the  energy  is 
so  small  relative  to  later  pulses,  these  synthetics  appear  to  have  delays. 

When  the  paths  of  the  potentials  become  shallower,  we  see  that  the  advance 
disappears.  The  planar  part  of  the  interface  begins  affecting  the  travel  times 
and  waveforms.  The  arrival  times  of  the  transmissions  through  the  bumps 
approach  those  of  transmissions  through  a  flat  interface.  The  exception  to  this 
pattern  is  the  residuals  of  the  plug  synthetics.  These  synthetics  have  an 
advance  of  .3  seconds  which  is  constant  as  a  function  of  horizontal  distance. 
This  behavior  results  from  the  constant  height  of  the  plug  across  the  entire 
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width  of  the  bump. 

What  is  the  relationship  between  the  travel  time  anomalies  and  the  ampli¬ 
tude  anomalies?  We  predict  that  as  the  magnitude  of  the  travel  time  anomalies 
decreases  the  magnitude  of  the  amplitude  anomalies  decreases  also.  The  ampli¬ 
tudes,  as  well  as  the  arrival  times,  will  be  controlled  by  the  planar  part  of  the 
surface.  This  relationship  is  observable  in  Figure  9.  The  amplitudes,  except  in 
one  case,  start  to  approach  the  value  of  .002  at  distances  ranging  from  5,000  to 
7,000  kilometers.  The  exception  is  the  maximum  and  first  pulse  amplitudes  of 
the  synthetics  from  the  plug  topography.  These  values  appear  to  systematically 
decrease  with  distance.  However  these  values  do  approach  the  planar  ampli¬ 
tudes  at  distances  beyond  7000  kilometers. 

Furthermore,  waves  which  arrive  earlier  than  is  predicted  by  planar  calcu¬ 
lations  also  have  lower  amplitudes  than  is  predicted.  Contrarily,  the  synthetics 
from  the  sine  topography  arrive  late  and  have  anomalously  high  amplitudes. 
Indeed,  the  pattern  of  residual  variation  is  precisely  mimicked  by  the  pattern  of 
amplitude  variation  for  this  structure.  The  mimicking  of  amplitude  and  travel 
time  anomalies  also  occurs  for  synthetics  from  the  exponential  bumps. 

This  mimicking  does  not  occur  for  synthetics  from  the  upwarp  and  plug 
topographies.  Each  of  these  synthetics  consists  of  multiple  arrivals.  Thus,  if  we 
take  the  maximum  amplitude  as  a  measure  of  amplitude  anomaly  and  compare 
with  travel  time  anomalies,  we  do  not  see  an  obvious  correlation  between  the  two 
parameters.  The  travel  time  is  perturbed  by  a  relatively  smell  part  of  the  sur¬ 
face.  The  maximum  amplitude  is  perturbed  by  a  much  larger  part  of  the  sur¬ 
face.  It  is  a  less  local  property  of  the  topography.  If  a  broader  band  source 
tune  function  interacts  with  the  surface,  the  amplitude  anomaly  would  change 
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but  not  the  travel  time  anomaly.  To  improve  the  correlation,  we  measure  the 
amplitude  of  the  first  pulse  of  the  synthetic  if  it  is  different  from  the  maximum 
amplitude.  These  values  are  shown  in  Figure  9  by  the  open  circles  and  triangles 
for  the  upwarp  and  plug  synthetics.  We  do  not  improve  the  visual  correlation 
significantly. 

The  modeling  of  a  symmetric  structure  demonstrates  that  such  a  structure 
on  the  koho,  consistent  with  travel  time  residuals,  causes  a  factor  of  2J£  in  max¬ 
imum  amplitude  variation  as  a  function  of  distance.  The  variation  of  amplitudes 
of  first  pulses  is  somewhat  less.  Neither  variation  is  as  large  as  the  observations 
of  amplitude  changes  between  stations  for  a  given  test  at  NTS.  Furthermore, 
where  there  is  a  large  amplitude  variation,  there  is  significant  waveform  distor 
tion  of  the  synthetics.  Also,  the  bumps  generally  cause  low  amplitudes.  We  note 
that  the  mesa  data  set  has  both  anomalously  high  and  low  amplitudes.  The  low 
amplitude  synthetics  arrive  early.  However,  any  relationship  between  amplitude 
anomalies  and  travel  time  anomalies  is  dependent  on  frequency  because  the 
effect  of  a  structure  on  a  wave  is  dependent  on  frequency. 

The  modeling  to  this  point  produces  anomalies  which  are  dependent  only  on 
distance  because  the  structures  are  symmetric.  We  now  introduce  asymmetry 
into  the  problem  by  allowing  the  source  to  move  off  the  center  of  a  symmetric 
upwarp.  We  do  these  calculations  because  observed  travel  time  and  amplitude 
anomalies  are  presented  as  a  function  of  azimuth.  Yet  we  do  not  know  if  these 
anomalies  arise  from  azimuthal  or  delta  heterogeneities.  Additionally,  there  is 
an  observed  variation  of  amplitude  at  a  given  station  writh  a  change  in  source 
position  in  the  mesa.  If  we  change  the  source  position  across  a  sample  struc¬ 
ture,  can  we  reproduce  the  factor  of  2%  seen  in  Figure  3?  We  also  wish  to 


examine  whether  there  is  any  systematic  relationship  between  amplitude  and 
travel  time  anomalies  as  a  function  of  azimuth  as  we  have  done  previously  with 
these  parameters  as  a  function  of  distance. 

The  modeling  experiment  is  similar  to  the  previous  one.  The  receivers  are 
1000  to  7000  kilometers  horizontally  away  from  the  center  of  the  topography 
and  20,000  kilometers  below  the  source.  The  sources  are  45  kilometers  above 
the  baseline  of  the  Moho.  To  produce  the  azimuthal  anomalies  in  the  synthetics, 
we  move  the  source  off  the  center  of  the  upwarp  in  one  direction  in  increments 
of  two  kilometers.  The  responses  are  calculated  for  seven  distances  at  five 
different  azimuths.  We  select  the  topography  upwarp  for  this  experiment.  The 
choice  of  topography  is  somewhat  arbitrary;  however,  we  make  this  particular 
choice  because  this  topography  causes  substantial  variation  in  amplitudes  as  a 
function  of  distance.  If  this  topography  fails  to  produce  much  azimuthal  varia¬ 
tion.  then  the  other  topographies  will  fail  to  do  so,  also. 

Figure  10  shows  our  results.  A  cross  section  of  the  source-receiver 
configuration  and  the  geometry  of  the  upwarp  is  in  the  center  of  the  figure.  In 
addition,  a  topographic  map  of  the  center  portion  of  the  grid  is  displayed.  The 
contours  are  in  kilometers  and  the  maximum  height  of  the  bump  is  10  kilome¬ 
ters.  The  topography  map  also  shows  the  source  locations  and  the  azimuthal 
lines  along  which  the  calculations  are  done.  The  resulting  waveforms  and  max¬ 
imum  amplitudes  surround  these  diagrams.  Each  group  of  28  waveforms  is  cal¬ 
culated  for  the  corresponding  azimuth.  The  groups  are  made  up  of  four  columns 
of  synthetics  corresponding  to  calculations  done  with  the  source  location  desig¬ 
nated  at  the  top  of  each  column.  Sources  A,  B,  C,  and  D  are,  respectively,  2,  4,  6 
and  8  kilometers  from  the  center.  The  rows  correspond  to  calculations  done  at 


the  horizontal  distances  next  to  the  row.  In  Figures  11,  12,  and  13  we  plot  the 
travel  time  residuals,  peak  amplitudes,  and  first  pulse  amplitudes  obtained  from 
these  synthetics  as  a  function  of  azimuth  for  each  distance. 

There  is  a  change  in  overall  complexity  of  the  synthetics  as  a  function  of 
azimuth.  The  waveforms  from  the  group  at  £?=0®  are  simple  and  impulsive  with 
relatively  high  amplitudes.  Only  the  stations  at  1,000  or  2,000  kilometers  have 
multiple  arrivals.  As  we  rotate  counter-clockwise  around  the  structure,  a 
greater  number  of  the  synthetics  in  each  group  have  multiple  arrivals,  and  con¬ 
sequently,  low  amplitudes.  Synthetics  at  6-  135°  and  180°  all  have  multiple 
arrivals.  The  reason  for  this  trend  is  the  same  as  in  the  previous  modeling  study. 
As  the  sources  move  in  the  direction  of  a  line  along  6-  0°.  a  greater  proportion 
of  the  elements  which  constitute  the  planar  part  of  the  grid  contribute  to  the 
potentials  calculated  in  the  direction  of  this  line.  Hence,  synthetics  of  this  line 
become  more  impulsive  as  the  source  migrates  from  position  A  to  position  D.  In 
contrast  the  synthetics  at  8 =  135®  and  180®  remain  complex.  The  elements 
which  contribute  to  these  potentials  are  largely  from  the  non-planar  part  of  the 
boundary. 

We  examine  the  maximum  amplitudes,  first  pulse  amplitudes,  and  travel 
time  residuals  in  Figures  11,  12,  and  13  for  systematics  as  a  function  of  distance 
or  azimuth.  The  behavior  of  maximum  amplitude  with  distance  and  azimuth  is 
the  most  variable  of  the  three  parameters.  The  maximum  amplitudes  as  a  func¬ 
tion  of  azimuth  do  not  correlate  very  well  with  the  travel  time  residuals.  The 
rapid  change  of  this  parameter  with  azimuth  and  source  position  reflects  the 
sensitivity  of  maximum  amplitudes  to  slight  changes  in  relative  timing  between 
the  two  arrivals  which  make  up  the  synthetics.  The  least  variation  of  maximum 
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amplitudes  with  distance  and  azimuth  occurs  at  synthetics  calculated  with 
source  position  A,  the  closest  source  to  the  center  of  the  symmetric  source.  We 
increase  this  variation  with  azimuth  and  distance  when  we  remove  the  source 
further  away  from  the  center  to  positions  B.  C.  and  D. 

The  maximum  amplitudes  at  6  =  135°  for  source  positions  B,  C.  and  D  are 
lower  than  the  amplitudes  at  other  azimuths.  This  trend  is  a  result  of  1)  the 
degree  to  which  the  planar  part  of  the  grid  contributes  to  the  response  and  2) 
the  degree  of  symmetry  of  the  source  and  receiver  locations  with  respect  to  the 
structure.  Stations  along  0=0°  and  8-  180“  are  in  positions  of  symmetry  with 
respect  to  the  sources.  Elements  on  either  side  of  a  line  which  divides  the  grid 
contribute  simultaneously  to  the  responses  at  these  stations  and  .consequently, 
cause  higher  maximum  amplitudes.  Stations  along  8-  90°  and  8-  45°  are  not 
symmetrically  positioned  with  respect  to  the  source;  however,  the  planar  part  of 
boundary  largely  contributes  to  these  responses;  thus,  they  have  high  maximum 
amplitudes.  But  receivers  along  6  =  135'  are  placed  asymmetrically  which 
causes  elements  to  illuminate  at  different  times;  in  addition,  these  elements  are 
largely  in  the  perturbed  part  of  the  boundary.  These  two  factors  combine  to 
produce  the  overall  lower  maximum  amplitudes  of  receivers  at  8=  135°. 

Although  we  do  not  discern  any  relationship  between  the  maximum  ampli¬ 
tudes  and  the  travel  time  residuals  as  a  function  of  azimuth,  we  see  a  correla¬ 
tion  between  the  first  pulse  amplitudes  and  the  travel  time  residuals.  Synthet¬ 
ics  which  have  a  constant  first  amplitude  as  a  function  of  azimuth  also  have 
approximately  constant  travel  Lime  residuals.  When  the  travel  time  advances 
increase  as  a  function  of  azimuth,  the  first  amplitudes  decrease  with  azimuth. 
Thus,  early  synthetics  have  lower  first  amplitudes  than  the  later  synthetics.  The 


trend  of  early  arrivals  with  low  amplitudes  and  late  arrivals  with  high  amplitudes 
holds  true  for  all  azimuths,  distances,  and  source  positions. 

The  travel  time  residuals  decrease  as  a  function  of  distance  at  all  azimuths 
except  for  8  -  135°  and  8  —  180°  for  sources  C  and  D.  At  these  azimuths,  the 
residuals  increase  as  the  distance  increases.  Clearly  if  we  pull  the  source  off  the 
center  far  enough,  the  shallower  rays  will  interact  with  the  upwarped  part  of  the 
topography  while  the  steeper  rays  interact  with  the  flat  part  of  the  grid. 

As  the  source  moves  off  the  center,  the  range  of  variation  of  first  ampli¬ 
tudes  and  travel  time  residuals  as  a  function  of  azimuth  exceeds  the  range  of 
these  parameters  as  a  function  of  distance;  that  is.  the  trends  of  these  parame¬ 
ters  are  stronger  in  azimuth  than  in  distance.  Thus,  stronger  variation  of  travel 
time  and  amplitude  anomaly  with  azimuth  than  with  distance  may  be  an  indica¬ 
tor  of  lateral  variation  with  azimuth,  as  well  as  with  distance,  despite  an  uneven 
station  distribution  of  the  existing  data  set. 

What  conclusions  can  be  drawn  from  this  modeling  experiment?  Firstly,  we 
create  a  variation  of  of  first  amplitudes  as  a  function  of  azimuth  and  source 
position.  The  change  of  amplitude  with  source  position  is  largest  at 
8-  0°  and  45°  and  is  the  least  at  8=  180®.  However,  the  variation  at  0®  and  45® 
is  created  at  the  cost  of  considerable  distortion  of  the  waveform.  This  feat  ure  of 
low  amplitude  waveforms  with  complex  or  broadened  pulses  and  high  amplitudes 
waveforms  with  simple  narrow  pulses  is  not  apparent  in  the  mesa  data  set. 

Secondly,  we  create  a  trend  in  the  first  pulse  amplitudes  with  azimuth. 
Specifically  we  cause  high  amplitudes  at  0=0®,  135®  ,  and  180®.  If  the  source 
moves  far  enough  away  from  upwarp,  it  causes  no  amplitude  anomalies.  How¬ 
ever  structure  can  produce  a  systematic  azimuthal  trend  in  amplitudes.  But  we 
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must  be  cautious  about  pushing  this  interpretation  too  far.  The  azimuthal  trend 
is  an  artifact  of  the  moveout  of  the  source  in  one  direction  with  respect  to  the 
lines  of  receivers.  If  we  were  to  distribute  sources  all  over  the  upwarp  and  then 
calculate  the  averages  of  the  first  amplitudes  at  each  azimuth  for  all  the 
sources,  we  would  undoubtedly  eliminate  any  trend  with  azimuth.  Thus,  the  sta¬ 
bility  of  the  amplitude  pattern  of  all  mesa  events  with  azimuth  location  is  not 
easily  explained  by  structure  on  the  Moho  or  any  unusual  velocity  plug  unless 
the  sources  are  fortuitously  located  to  one  side  of  the  heterogeneity. 

Thirdly,  we  see  a  visual  correlation  between  travel  time  residuals  and  ampli¬ 
tudes  of  first  pulses,  but  do  not  see  any  between  residuals  and  peak  amplitudes. 
This  correlation  may  be  diagnostic  of  structure  as  opposed  to  tectonic  release. 

Discussion  and  Conclusions 

The  previous  two  modeling  experiments  show  that  a  structure  on  the  Moho 
which  causes  travel  time  residuals  compatible  with  the  Minster  et  al.,  (1981)  and 
the  Spence  (1974)  studies  can  produce  variations  of  2}£  of  amplitudes  as  a  func¬ 
tion  of  delta,  azimuth,  and  source  position.  The  variation  is  created  at  the  cost 
of  considerable  distortion  of  the  waveform.  Furthermore  the  travel  time  residu¬ 
als  correlate  with  first  pulse  amplitudes  but  not  with  peak  amplitudes.  To  see 
whether  these  initial  results  are  relevant  to  the  Pahute  Mesa  waveform  data,  we 
must  now  include  a  pP  phase,  a  Q  operator  and  a  short  period  instrument  opera¬ 
tor  in  a  few  Kirchhoff  synthetics. 

We  choose  two  sets  of  five  KirchhofT-Helmholtz  synthetics  calculated  at  4000 
kilometers  previously  for  the  azimuthal  study.  The  two  sets  correspond  to  the 
source  positions  A  and  D  at  five  azimuths.  This  choice  represents  two  extremes 


of  source  positions  relative  to  a  structure  and  may  give  us  a  reasonable  idea  of 
what  to  expect  in  amplitude  and  waveform  variation  as  test  sites  move  within  the 
mesa 

We  put  pP  into  the  Kirchhoff-Helmholtz  synthetics  by  convolving  these  syn¬ 
thetics  with  a  boxcar  function  of  unit  height  and  a  width  corresponding  to  the 
pP-P  lag  time.  This  convolution  yields  the  impulse  response  of  P  and  pP  if  the 
incident  source  is  a  modified  Haskell  function  rather  than  its  time  derivative. 
We  justify  this  simple  model  of  pP  by  assuming  that  the  reflection  coefficient  of 
this  phase  is  -1.  Lay  et  al.  (1933a)  estimates  the  reflection  coefficient  as  .96  for 
pP.  We  further  argue  that  pP  interacts  with  the  same  part  of  the  surface  as  does 
P.  This  assumption  is  good  for  the  shallow  depths  of  the  rnesa  tests  which  rangr 
from  .5  to  1.4  kilometers.  The  width  of  the  boxcar  is  .85  seconds;  this  estimate  of 
the  pP  -  P  lag  time  is  taken  from  Lay  et  al.  (1983).  A  short  period  instrument  and 
a  Futterman  Q  operator  is  also  convolved  into  these  synthetics.  WTe  use  a  Haskell 
function  with  parameters  (D  =  2,K  =  10)  while  Lay  et  al.  (1983a)  use  slightly 
different  values  (B  =  1,K  =  8). 

Figure  14  displays  the  results  of  the  convolutions.  The  first  column  con¬ 
tains  the  initial  Kirchhoff-Helmholtz  synthetics  with  peak  amplitudes  taken  from 
Figure  10.  The  second  column  show  these  synthetics  convolved  with  a  WWSSN 
short  period  instrument  and  a  boxcar.  We  introduce  additional  complexity  into 
the  waveform  but  do  not  change  the  range  of  peak  amplitudes  significantly  when 
we  include  an  instrument  and  pP.  The  complexity  of  the  waveform  caused  by 
structure  is  masked  by  the  dominant  interference  between  P  and  pP. 

We  next  convolve  these  synthetics  with  Futterman  Q  operators  with  a  t*  of 
.5  and  1.  The  waveforms  and  their  ab  amplitudes  are  displayed  in  columns  3  and 


4.  The  ab  amplitudes  are  plotted  as  function  of  azimuth  for  both  sources  and  t* 
values  in  Figure  15.  We  remove  the  complexity  of  the  waveform  for  both  sources 
with  the  two  t*  values.  However  there  are  some  observable  differences  in  the 
first  and  third  peaks  of  the  waveforms  as  a  function  of  azimuth.  The  first  peak 
widens  as  azimuth  increases.  The  third  peak  becomes  smaller  and  disappears 
altogether.  Moreover.  Figure  15  shows  a  variation  of  ab  amplitudes  with  azimuth 
of  2 for  source  D  if  t*  is  .5.  However,  when  t*  is  1.  this  variation  reduces  to  a 
factor  of  1.7.  We  also  obtain  a  variation  of  2$£  of  ab  amplitudes  with  respect  to 
source  position  if  t*  is  .5.  This  occurs  at  0=0°  and  0=45°;  however  the 
difference  between  ab  amplitudes  for  the  two  sources  decreases  as  the  azimuth 
increases. 

Thus  we  cannot  predict  the  observed  ab  amplitude  variation  with  azimuth 
or  station  location  if  we  use  a  structure  on  the  Moho  10  kilometers  high  and 
approximately  25  kilometers  wide.  If  the  factors  of  12  or  40.  seen  in  Figures  3, 
4,  and  5,  are  measures  of  a  purely  near-source  phenomena,  then  one  requires  a 
structure  several  hundred  kilometers  in  scale  on  a  boundary  to  match  these  fac¬ 
tors.  This  structure  would  distort  the  waveform  considerably.  Yet  there  is  no 
obvious  evidence  in  the  observed  seismograms  for  a  correlation  between  low 
amplitudes  and  complicated  and/or  broadened  waveforms  or  high  amplitudes 
and  simple,  impulsive  waveforms.  We  speculate  that,  rather  than  a  large  struc¬ 
ture  on  a  single  boundary,  a  small  velocity  or  density  perturbation  along  a 
several  hundred  kilometer  ray  path  may  produce  the  desired  amplitude  change. 
However  we  cannot  test  this  speculation  with  our  method. 

On  the  other  hand,  we  predict  a  factor  of  2%  in  ab  variation  with  source 
location  if  t*  is  .5  .  This  variation  is  not  accompanied  by  any  significant 


y 


v-  -.iv- "--  *  •• 

'& 


±2- 


112 


waveform  distortion.  Although  we  only  produce  a  factor  of  2^  at  two  azimuths, 
this  is  an  artifact  of  the  source  moveout  across  the  structure.  If  sources  were 
uniformly  distributed  over  the  structure,  we  would  cause  this  same  magnitude 
of  variation  at  all  the  azimuths.  Furthermore,  no  source  would  be  systematically 
higher  in  ab  amplitudes  than  another  source  at  all  the  azimuths.  Unfortunately 
Figure  3  does  not  show  systematics  with  respect  to  source  location.  The  data 
should  be  examined  for  such  trends. 


Also  the  travel  time  residuals  do  not  correlate  with  the  ab  patterns.  The 
convolution  with  two  t*  values  demonstrates  the  frequency  dependence  of  the 
phenomena.  If  t*  is  larger  than  1,  we  will  produce  a  flat  pattern  of  ab  amplitudes 
with  azimuth  but  the  travel  time  residuals  will  not  change.  Thus  we  do  nol 
expect  a  systematic  relationship  between  travel  time  and  amplitude  anomalies 
because  the  ab  amplitudes  are  sensitive  to  Q  but  the  travel  times  are  relatively 
stable.  This  is  an  unfortunate  result  because  such  a  correlation  would  be  diag¬ 
nostic  of  structure  as  opposed  to  tectonic  release.  As  yet,  there  has  been  no 
study  which  definitively  demonstrates  a  relationship.  In  addition  the  data  sets 
of  amplitude  and  travel  time  measurements  do  not  have  a  one-to-one  correspon¬ 
dence.  Lay  et  al.  (1983a)  measure  the  ab  amplitudes  off  of  short  period  YfWSSN 
instruments  while  Minster  et  al.  (1981)  and  Spence  (1974)  use  culled  travel  time 
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measurements  from  the  ISC  catalogs.  We  clearly  need  a  study  which  compares 
the  travel  time  and  amplitude  from  the  same  seismogram. 

The  largest  overall  variation  in  ab  amplitude  with  azimuth  in  Figure  15 
occurs  because  of  differences  in  t*  values,  Tet  there  are  no  noticeable  changes 
in  the  waveform.  Thus  the  ab  amplitudes  are  far  more  sensitive  than  waveforms 
to  Q  differences.  Perhaps,  a  lateral  variation  of  Q  with  path  can  produce  the 
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extreme  scatter  of  ab  amplitudes  for  both  tests  inside  and  outside  the  mesa. 
However  it  cannot  explain  the  differences  in  patterns  between  these  different 
test  site  areas.  If  the  near-mesa  anomalies  are  .indeed,  at  least  a  factor  of  10, 
then  structure  on  the  Moho  which  is  compatible  with  travel  time  residuals  can¬ 
not  produce  these  large  variations  of  amplitude  with  azimuth.  However  such  a 
structure  could  explain  the  observed  variation  of  ab  amplitudes  with  source 
position  at  a  given  station. 
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Figure  Captions 


la.  The  geometry  of  the  Kirchhoff-Helmholtz  calculations  for  transmission 

across  two  acoustic  media  with  sound  speeds  Oj  and  a2  and  densi¬ 
ties  Pi  and  p2.  The  source  is  in  Vx  at  Xj,  and  the  receiver  is  in  F2  at 
JL- 

lb.  A  close  up  of  a  piece  of  the  boundary  which  displays  the  angles. 

2a.  Two  synthetics  and  the  grid  geometry  used  to  compute  them.  Synthetic  A  is 
contaminated  by  a  truncation  phase  which  originates  from  the 
edge  of  the  grid.  Synthetic  B  is  contaminated  by  a  phase  which 
originates  from  the  abrupt  change  in  boundary  conditions.  The 
grid  next  to  synthetic  B  is  gray  when  tp  -  0  on  the  boundary. 

2b.  A  comparison  between  Kirchhoff-Helmholtz  and  first  motion  solutions.  The 
input  source  is  the  first  derivative  of  a  modified  Kaskell  function 
with  parameters  (B=2,K=10).  The  maximum  dimensionless  ampli¬ 
tude  of  the  source  input  function  is  45.1  . 


3.  The  short  period  P  wave  ab  amplitude  data  set  for  25  Pahute  Mesa  events 
plotted  as  a  function  of  station  location.  The  amplitudes  are 
corrected  for  event  size,  geometric  spreading  and  instrument  gain 
at  1  second  and  are  plotted  relative  to  a  master  event  (from  Lay  et 
al..  1983a). 


4.  The  relative  ab  amplitudes  of  GREELEY  and  FAULTLESS  as  a  function  of  sta¬ 

tion  location  (from  Lay  et  al..  1983a). 

5.  Ratios  of  relative  ab  amplitudes  of  FAULTLESS.  PILEDR1VER.  and  B1LBY 

divided  by  the  average  relative  ab  amplitudes  of  the  mesa  events 
(from  Lay  et  al..  1983b). 

6.  Seismograms  from  FAULTLESS  (top  record)  and  GREELEY  (bottom  record) 

displayed  in  order  of  increasing  azimuth  in  the  range  of  30°  to 
100°.  Also  shown  are  the  absolute  ab  amplitudes  in  millimicrons, 
corrected  only  for  instrument  gain  at  1  second. 

7.  Transmitted  potentials  from  sources  35  kilometers  above  the  center  of  the 

structure.  The  cross  sections  of  the  structures  upwarp,  exponen¬ 
tial,  plug  and  sine  are  above  the  synthetics.  For  comparison 
potentials  which  propagate  through  a  flat  boundary  are  shown  in 
the  first  column.  The  potentials  are  from  receivers  which  are 
20,000  kilometers  below  the  source  and  which  vary  from  0  to  7000 
kilometers  horizontally  away  from  the  source.  All  amplitudes  are 
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multiplied  by  .01. 


8.  Travel  time  contours  for  a  source  35  kilometers  directly  above  the  structure 
and  receivers  20,000  kilometers  directly  below  the  source.  The 
four  structures  are  a)  a  plane  ,  b)  an  upwarp,  c)  a  plug,  and  d)  a 
sine  function.  The  contours  are  projected  onto  the  topographies 
and  flat  grids.  The  synthetics  which  correspond  to  each  travel 
time  projection  are  in  between  the  two  projections.  The  contour 
interval  is  .125  seconds  as  are  the  tick  marks  of  the  synthetics. 
The  geometric  arrival  time  is  the  center  of  the  contours. 


9.  Plots  of  peak  amplitudes,  amplitudes  of  first  pulses,  and  travel  time  residuals 
as  a  function  of  distance  from  synthetics  in  Figure  7.  The  different 
symbols  correspond  to  different  topographies  and  are  at  the  bot¬ 
tom  of  the  figure.  Where  first  pulse  amplitudes  are  different  from 
peak  amplitudes,  the  values  of  first  pulse  amplitudes  are  plotted 
with  open  symbols  and  the  peak  amplitudes  are  plotted  with  closed 
symbols.  The  dotted  line  corresponds  to  the  peak  amplitudes  from 
synthetics  which  propagate  through  a  planar  interface. 


10.  Synthetics  from  the  topography  upwarp  calculated  for  four  source  positions, 
five  azimuths  and  seven  distances.  The  topography  map  and  cross 
section  with  source  positions  are  in  the  center.  The  contour  inter¬ 
val  is  1  kilometer.  The  distances,  angles,  and  azimuths  of  the 
receivers  are  also  shown. 


11.  Travel  time  residuals  for  source  locations  A,  B.  C,  and  D  plotted  as  a  function 
of  azimuth  and  distance. 


12.  Peak  amplitudes  for  source  locations  A.  B.  C,  and  D  plotted  as  a  function  of 
azimuth  and  distance. 


13.  Amplitudes  of  the  first  pulse  for  source  locations  A,  B,  C,  and  D  plotted  as  a 
function  of  azimuth  and  distance. 


14.  Kirchhoff-Helmholtz  responses  (first  column)  convolved  with  a  boxcar  of 
width  .85  seconds  and  a  short  period  WWSSN  instrument  (second 
column)  and  a  Q  operator  with  t*  values  of  .5  (third  column)  and  1 
(fourth  column).  Responses  are  from  a  distance  of  4000  kilome¬ 
ters,  azimuth  range  of  0®  to  180®  and  source  locations  A  and  D. 


15.  ab  amplitudes  from  synthetics  in  Figure  14  plotted  as  a  function  of  azimuth 
for  both  source  locations  A  and  D  and  bot  h  t®  values. 
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are  required  to  construct  synthetics  for  mixed  path  situations.  This  is 
accomplished  by  allowing  locally  dipping  structure  and  making  some 
modifications  to  generalized  ray  theory.  ]  Local  ray  parameters  are  expressed 
in  terms  of  a  global  reference  which  allows  a  de-Hoop  contour  to  be 
constructed  for  each  generalized  ray  with  the  usual  application  of  the 
Cagniard  de-Hoop  technique.  Several  useful  approximations  of  ray  expansions 
and  WKBJ  theory  are  presented.  Comparisons  of  the  synthetics  produced  by 
these  two  basic  techniques  with  known  solutions  demonstrates  their  reliability 
and  limitations. 

In  Section  III,  we  have  modeled  the  SH  motion  from  earthquakes  in  the 
northwest  Atlantic  ocean  to  derive  an  oceanic  upper  mantle  shear  velocity 
model.  The  signals  were  recorded  on  long-period  WWSSN  and  Canadian  network 
stations'  on  the  east  coast  of  North  America.  This  data  indicates  a  fast  (4.75 
km/sec)  lid  of  about  100  km  thickness  in  the  older  western  Atlantic.  Given 
the  lid  structure,  the  waveforms  and  traveltimes  from  the  more  distant  data 
put  tight  constraints  on  the  shear  velocities  at  greater  depths.  The  velocity 
below  20Q  km  was  found  to  be  indistinguishable  from  a  model  of  the  East 
Pacific  Rise  (Grand  and  Helmberger,  1983)  found  using  the  same  technique.  We 
find  the  Canadian  shield  to  be  faster  than  both  the  old  northwestern  Atlantic 
and  the  young  East  Pacific  Rise  to  about  400  km  depth.  No  variations  below 
400  km  are  necessary  to  explain  the  data. 

•  In  Section  IV,  we  extend  the  Kirchhof f-Helmholtz  integral  method  to 
calculate  acoustic  potentials  which  transmit  through  three  dimensional  warped 
boundaries.  *  We  specify  the  potentials  on  an  arbitrary  surface  with  Snell's 
law  and  plane-wave  transmission  coefficients  and  numerically  integrate  their 
contributions  at  a  receiver  via  the  scalar  integral  representation  theorem. 

The  method  is  appropriate  for  modeling  precritical  transmitted  potentials. 
Results  from  test  models  compare  well  with  optical  solutions  for  transmissions 
through  a  flat  interface. 


